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ABSTRACT 


The  Naval  Post  Graduate  School  has  obtained  software  which  can  be  used  for  the 
thermal  analysis  of  an  electronic  component.  This  software  includes  two  versions,  one  for 
steady  state  analysis  and  the  other  for  transient  analysis.  Each  version  consists  of  two 
programs,  a  model  builder  and  a  thermal  analyzer.  This  thesis  describes  a  user  friendly, 
menu  driven  specific  model  builder  which  may  be  used  to  rapidly  generate  a  thermal  model, 
containing  up  to  750  nodes,  for  a  microcircuit  die.  This  model  builder  is  an  improvement 
over  the  existing  model  builder  in  that  the  only  pertinent  input  will  be  the  physical 
dimensions  and  heat  transfer  data.  The  main  feature  of  this  model  builder  is  the 
accomodation  of  the  maiked  variation  of  silicon  thermal  conductivity  with  temperature. 
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I.  INTRODUCTION 


Microcircuit  chips  have  become  commonplace  items  in  all  types  of  electronic 
equipment,  both  military  and  commercial.  The  capabilities  of  these  chips  continue  to 
increase  as  improvements  in  design  and  manufacturing  are  incorporated.  However,  the 
design  process  becomes  more  difficult  as  the  ability  to  place  more  and  more  circuits  on 
these  chips  is  realized.  One  area  of  immediate  concern  is  in  thermal  control.  Chips 
containing  more  than  one  million  components  in  a  60  mm^  area  are  known  to  exist  [Ref. 
1].  Heat  dissipation  from  components  contained  in  the  chip  as  well  as  from  any  outside 
heat  sources  needs  to  be  thoroughly  analyzed.  The  reasons  for  conducting  a  thermal 
analysis  on  a  microcircuit  structure  concern  the  impact  of  high  temperature  on  reliability, 
bias  stabilization,  manufacturing  costs,  and  thermal  runaway. 

The  reliability  of  a  microcircuit  is  directly  related  to  its  junction  temperature.  It  is  now 
common  knowledge  that  failure  rate  is  proportional  to  junction  temperature  and  designers 
plan  for  the  microcircuit  to  operate  within  a  specifled  temperature  range.  If  this  range  is 
exceeded,  the  operating  capability  of  the  device  may  be  severely  reduced  or  it  may  fail 
completely.  Because  the  failure  rate  increases  exponentially  as  the  temperature  increases, 
the  effects  of  high  junction  temperature  can  be  drastic. 

It  is  necessary  to  bias  an  electronic  circuit  in  order  to  obtain  a  stable  and  predictable 
operating  point.  As  the  junction  temperature  changes,  the  impact  on  this  operating  point 
will  be  directly  impacted.  Consider  Figure  1  and  observe  that  by  Kirchhoffs  voltage  law, 

-Vcc  +  Rclc  +  Vce  +  R«lc  =  0  (1) 

so  that  the  collector  to  emitter  voltage  will  be 

Vce  =  Vcc-Ic(Rc  +  Re)  (2) 
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Vcc  =  power  supply  voltage  to  collector 
B  =  base 
C  =  collector 
E  =  emitter 

Vce  =  voltage  across  collector  and  emitter 
Rg  =  emitter  resistance 
Rc  =  collector  resistance 
le  =  emitter  current 
I  ^  =  collector  current 


It  is  a  fact  that,  the  collector  current,  Ic,  is  directly  proportional  to  the  junction  temperature. 
As  the  junction  temperature  increases,  the  voltage  drop  across  Rc  and  increases  and  as  a 
result,  Voe  is  quite  likely  to  decrease  to  an  intolo^ly  low  level. 

During  manufacturing,  the  chip,  at  some  point,  must  undergo  a  soldering  process  in 
which  connections  are  made  from  chip  to  chip  or  from  one  component  to  another.  A 
particular  soldering  compound  must  be  selected  which  will  accomodate  the  expected 
operating  temperature  range  without  failure.  The  soldering  compound  cost  is  based  on  it’s 
composition  and,  as  general  rule  of  thumb,  the  higher  the  melting  point,  the  higher  the  cost 
An  increased  cost  incurred  by  a  manufacturer  for  a  high  temperature  solder  may  be  very 
minimal  on  a  per  unit  basis.  But  any  manufacturer  dealing  with  microcircuit  chips  may  be 
handling  millions  of  them  and  these  costs  could  conceivably  reach  the  million  dollar  mark 
which  is  probably  not  recoverable. 

Hot  spots  anywhere  in  the  microcircuit  chip  could  result  in  thermal  runaway.  Because 
the  temperature  increases  as  the  collector  current  increases,  a  hot  spot  near  the  transistor 
will  cause  an  increase  in  the  collector  current.  This  then  causes  an  increase  in  power 
dissipation  (P  =  IcVc*)  which  in  turn,  causes  an  increase  in  temperature.  This  increase  in 
temperature  causes  an  additional  increment  of  collector  current  to  flow  and  leads  to  yet 
another  increase  in  power  dissipation.  This  regenerative  process  is  known  as  thermal 
runaway  and  as  the  cycle  continues  it  could  lead  to  drastic  and  complete  device  destruction. 
[Ref.  1] 

The  Naval  Post  Graduate  School  has  obtained  software  which  can  be  used  for  the 
thermal  analysis  of  an  electronic  component.  This  software  includes  two  versions,  one  for 
steady  state  analysis  and  the  other  for  transient  analysis.  Each  version  consists  of  two 
programs,  a  model  builder  and  a  thermal  analyzer.  The  model  builder  is  designed  to 
provide  the  user  with  the  means  to  model  a  physical  configuration  in  order  to  provide  the 
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thermal  analyzer  with  an  input  file.  This  file  is  taken  by  the  thermal  analyzer  as  its  input 
and  an  output  file  containing  a  summary  of  the  temperatures  within  the  configuration  will 
be  generated. 

The  generation  of  the  thermal  model  is  the  biggest  hurdle  in  the  actual  analysis  and  the 
analyst  is  powerless  in  the  absence  of  a  model  builder  program.  However,  the  model 
building  process  begins  with  a  drawing  of  the  configuration  and  a  subdivision  of  it  into 
small  but  finite  subvolumes.  Each  subvolume  is  presumed  isothermal  and  the  centers  of 
each  subvolume  are  then  representative  of  the  entire  subvolume.  These  centers  (centroids) 
are  called  node  points  or  merely  nodes  and  thermal  analysis  involves  the  writing  of  n  node 
equations  in  n  unknown  temperatures  where  the  nodes  are  connected  by  thermal 
conductances. 

The  very  laborious  and  time  consuming  process  of  writing  the  node  equations  can  be 
avoided  through  the  use  of  a  model  builder  which  will  write  them  automatically.  Yet,  there 
are  drawbacks  in  the  use  of  a  general  model  builder  in  that,  while  it  can  be  used  for  a 
specific  geometry  and  set  of  material  conditions,  many  hours  can  be  spent  on  input  in  front 
of  a  computer  terminal  and,  before  using  any  model  builder,  careful  consideration  must  be 
given  to  whether  or  not  the  configuration  under  analysis  will  be  sufficiently  well 
represented  by  the  model  generated. 

This  thesis  describes  a  specific  model  builder  which  may  be  used  to  rapidly  generate  a 
thermal  model  containing  up  to  750  nodes  for  a  microcircuit  die.  The  main  feature  of  this 
model  builder  is  the  accomodation  of  the  marked  variation  of  silicon  thermal  conductivity 
with  temperature. 

The  specifrc  model  builder  has  additional  features.  These  will  be  described  in  this 
thesis  and  are  as  follows: 
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1)  The  capability  of  working  in  SI  or  English  units. 

2)  The  user  has  five  aspect  ratios  from  which  to  choose. 

3)  The  provision  of  up  to  750  nodes  fen*  the  die  model. 

4)  The  ability  to  input  heat  dissipatiem  using  several  methods. 

5)  The  automatic  calculation  of  conductance  values  based  on  user  input 
dimensions. 

6)  The  automatic  correction  for  variations  in  the  thermal  conductivity  values  of 
silicon  with  temperature  and  location. 

7)  The  provision  of  automatic  linearization  of  radiation. 
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II.  A  PROBLEM  IN  THERMAL  ANALYSIS 


The  provision  of  a  solution  to  a  simple  problem  in  heat  transfer  analysis  will  show  that 
a  thermal  analyzer  is  required  for  serious,  large  scale,  thermal  analysis. 

Conduction  is  the  heat  flow  mechanism  whereby  heat  is  transfered  from  one  part  of  a 
homogeneous  medium  under  the  influence  of  a  temperature  gradient  without  appreciable 
displacement  of  the  particles  that  compose  the  medium.  It  was  Fourier  who  proposed  that 
the  heat  flow  is  directly  proportional  to  the  area  of  the  heat  flow  path  and  the  temperature 


(3) 

(4) 


A  =  area  of  the  heat  flow  path 

k  =  thermal  conductivity  of  the  material 
cfT 

^  =  change  in  temperature  per  unit  path  length  (temperature  gradient) 
where  the  minus  sign  assures  a  positive  heat  flow  in  the  presence  of  the  required  negative 
temperature  gradient.  [Ref.  2] 

Equation  (4)  serves  to  define  the  thermal  conductivity  of  the  material 

k  =  -  W/m-oK  (5) 

^dx 


and  it  must  be  noted  that  the  simplicity  of  this  development  quickly  changes  when  kisa 
function  of  temperature,  as  in  the  case  of  silicon. 

Equation  (4)  is  a  simple  differential  equation  which  can  be  solved  using  the 
boundary  conditions  in  Figure  2.  Separation  of  the  variables  gives: 
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-idx=dT 

and  the  integration  between  the  limits: 

Xj  1 1 

gives: 

-j5c<*2-Xl)=T2-T, 

With  L  =  X2  -  xi,  it  is  seen  that: 

q  =  !^(Ti-T2) 

which  is  called  Fourier’s  law. 


(6) 


Consider,  Ohm’s  law  which  governs  the  flow  of  current  in  an  electrical  resistor 
AV 


1  = 


(7) 


AT<=>AV 
L 


and  note  the  similiarity  between  the  forms  of  equations  (6)  and  (7).  An  analogy  may  be 
proposed  where  the  analogous  quantities  are: 

q  <=>  I  (8a) 

(8b) 

l^«R  (8c) 

In  particular,  equation  (8c),  defmes  the  thermal  resistance  for  the  conduction  process  and 
one  observes  that,  for  the  slab  in  Figure  2,  if  Ti  and  T2  are  given,  then  q  can  be  determined 
if  the  thermal  resistance  of  the  slab  is  known.  On  the  other  hand,  if  Ti  and  q  are  known, 
T2  can  be  calculated  using  R  = 


T2  =  Ti-Rq 

Convection  is  a  fluid  flow  process  and,  if  it  applies  in  heat  transfer,  it  is  the  transfer  of 
heat  from  or  to  a  conflning  surface  by  a  flowing  fluid.  The  fluid  flow  may  be  induced  by 
buoyancy  or  density  gradients  within  the  fluid  and  this  phenomenon  is  called  natural  or  free 
convection.  If  the  fluid  flow  is  induced  by  a  fan,  a  pump,  or  a  blower,  it  is  termed  forced 
convection.  [Ref.  2] 
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Heat  flow  by  convection  is  directly  proportional  to  the  temperature  difference  between 
the  confining  surface  and  fluid  and  the  surface  area  over  which  the  process  takes  place: 

q  «:  A  AT  (9) 

and  a  proportionality  factor  may  be  inserted: 

q  =  h  A  AT  (10) 

Equation  (10)  is  the  Newton  law  of  cooling  and  h  is  called  the  coefficient  of  heat  transfer 
defined  by: 

h=— 3— W/m2-oK  (11) 

A  AT 

Notice  that  the  concept  of  thermal  resistance  can  be  extended  to  the  convective  case 
where: 

(12) 

If  convection  exists  on  both  faces  of  the  slab  in  Figure  2,  the  configuration  of  Figure  3 
exists  and  the  electric  analog  is  shown  in  Figure  4.Here  the  total  resistance  is: 

■  hi  A  ^  kA  ^  h2A 


hi  k 


and  the  heat  flow  is  seen  to  be: 


(13) 


Thus  far,  no  difficulty  has  been  encountered  in  a  calculation  procedure.  The  picture 
changes  dramatically  when  radiation  is  considered. 

Radiation  is  electromagnetic  in  nature  and  is  the  means  for  heat  communication  in  a 
vacuum.  It  is  beyond  the  scope  of  this  thesis  to  derive  the  general  radiation  equation: 
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2 


Figure  3 
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q  =  oFAFeA(Ts4-T,4) 


(14) 


where 

a  is  the  Stefan-Boltzmann  constant  =  5.67(10-*)  W/m^-K^ 

Fa  is  a  shape  factor  that  accounts  for  the  arrangement  of  the  radiating  source  and 
the  absorbing  receiver 

Fe  is  the  emissivity  factor  that  accounts  for  the  ability  of  the  surfaces  to  absorb  or 
emit  radiation 

Ts  is  the  temperature  of  the  source,  Kelvin 
Tr  is  the  temperature  of  the  receiver,  Kelvin 
A  is  the  surface  area,  m^ 

In  order  to  put  the  radiation  on  a  q  =  K  AT  basis  where  K  is  the  thermal  conductance, 
consider  that; 

Ts4-T4  =  (T*2+T,2)(Ts2-T,2) 

and: 

T/ .  T  4  =  (T2  +  t,2)  (T,  +  T,)  (Ts  -  Tr) 

Equation  (14)  may  then  be  written  as: 


q  =  o  Fa  Fe  A  (Ts2  +  T,2)  (T*  +  T,)  (T*  -  T,) 

(15) 

and  a  radiative  heat  transfer  coefficient  may  be  defined: 

q  =  h,A(Ts-T,) 

(16) 

where 

hr  =  aFAFe(Ts2+T,2)Cr,+T,) 

(17) 

A  thermal  resistance  in  radiation  can  then  be  proposed  as: 

R  =  r-^ 
hr  A 

(18) 

and,  as  shown  in  Figure  5,  this  resistance  is  in  parallel  with  the  convective  resistance. 
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A  solution  procedure  for  the  thermal  network  in  Figure  5  must  be  iterative  because  of 
the  form  of  h,  given  by  equation  (17).  This  is  obviously  no  longer  a  simple  problem  and, 
at  this  point,  one  must  turn  to  the  computer  for  help  with  more  complicated  configurations. 


14 


III.  FINITE  DIFFERENCES  AND  NODAL  ANALYSIS 


It  is  the  purpose  of  this  section  to  provide  the  assurance  that  the  method  of  modelling 
will  yield  accurate  results  and  to  show  how  the  n  node  equations  in  the  n  unknown 
temperatures  are  formulated. 

Equation  (4)  is  but  a  one  dimensional  simplification  of  the  more  general  equation  of 
heat  conduction  [Ref.  3]: 


where 

T  =  temperature,  ®F  or  «€ 
t  =  time,  hrs. 

qj  =  rate  of  internal  heat  generation,  Btu/hr  or  Watts 
k  =  thermal  conductivity,  Btu/hr-ft-°F  or  Watts/m-®C 
p  =  specific  weight,  Ibm/ft^  or  kg/m^ 

X,  y,  z  =  coordinate  system  distances,  ft  or  m 
c  =  specific  heat,  Btu/lbrn-^F  or  Watt-hr/kg-°C 
We  note  that  if  we  wish  to  formulate  a  problem  that  deals  in  discrete  Ax,  Ay,  Az,  and 
AT’s,  we  must  be  able  to  accurately  approximate  the  first  and  second  derivatives. 

In  general  we  must  be  able  to  determine  the  slope  of  a  line  tangent  to  the  point  of 
interest  on  the  plot  of  the  temperature  as  a  function  of  x,  y,  or  z.  (Figure  6).  Using  x  as  the 
spacial  coordinate,  we  observe  that  the  slope  of  the  line  (m),  if  it  exists,  will  be: 

rr.-  ^  /(Xq+Ax)  -  /(Xq) 

"■-Ax^O 

The  term,  [/(Xq+Ax)  -  /(Xo)]/Ax,  is  called  the  difference  quotient  and  is  the  ratio  of  the 
change  in  the  value  of  the  function  at  Xq  and  Xo+Ax  to  the  change  in  x  .  The  limit  of  the 
difference  quotient  is  called  the  derivative  of  the  function  at  x,,: 
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Figure  6 
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lim  /(Xq+Ax)  -  /(Xq) 

J  "  Ax^O  Ax  ^  ^ 

In  examining  Figure  6  we  can  derive  the  derivative  for  this  particular  curve  at  the  point  Xo 

by  letting  /(x)  =  T(x)  and  noting  that: 

/(Xo+Ax)  =  T(Xo+Ax) 

/(Xo)  =  T(Xo) 

\  ^  T(Xo+Ax)  -  T(Xo) 

•'^*°^"Ax-»0  Ax  ^  ^ 

and  when  Ax  — »  0  we  observe  that  the  expression  for  the  first  derivative  becomes: 

/'(Xo)  (23) 

Ax 

and  we  note  that: 

/'(xo)=^  (24) 

Ax 

as  long  as  Ax  is  small  as  required  by  the  limiting  process.  [Ref.  4] 

The  second  derivative,  using  the  defintion  of  a  derivative,  will  be  the  ratio  of  the 
change  in  the  value  of  the  first  derivative  to  the  change  in  Ax.  Thus,  the  second  derivative 


in  Figure  6  becomes: 

P(  \  T(Xo+Ax)  -  T(Xo) 

Ax 


/'(Xo+Ax)  = 


T(Xo+2Ax)  -  T(Xo+Ax) 


„  lim  /'(Xq+Ax)  -  /'(Xq) 

J  ~  Ax-»0  Av 


r(xo)= 


T(Xo+2Ax)  -  T(xo+Ax)  T(xo+Ax)  -  T(xo) 


f '(Xo)  = 


T(xo+2Ax)  -  2T(Xo+Ax)  +  Tfxp) 
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As  in  the  case  of  the  first  derivative,  the  second  derivative  given  by  equation  (25)  is  a 
valid  approximation  as  long  as  Ax  is  small.  The  use  of  a  Taylor  series  expansion  will 
confirm  the  validity  of  equation  (25).  Taylor’s  theorem  states  that  a  function  can  be 
approximated  by  a  polynomial  of  the  form: 

/(X)  =  /(a)  +  (a)(x-a)  +  ^x-a)2  + ...  +  ^-^x-a)"  (26) 

where  the  polynomial  is  for  the  function  expanded  about  x  =  a  [Ref.  4].  Using  our 

example  the  Taylor  series  for  the  function  T(x)  at  T(Xo+2Ax)  is: 

rr./  A  N  rr,/  AX  dT(Xo+Ax)  .  1  d2T(xo+Ax)  . 

T(Xo+2Ax)  =  T(xo+Ax)  + - ^ — -  Ax  +  2 - ^[^^2 - 


1  d3T(xo+Ax) 


Ax^  +  ...  + 


1  d"T(Xo+Ax) 


6  dx^ 
and  for  the  function  T(x)  at  T(Xo): 


n! 


dx" 


Ax" 


X  .AX  dT(xo+Ax)  ,  .  1  d2T(Xo+Ax)  , 

T(Xo)  =  T(Xo+Ax)  - - 37 - Ax  +  TT - — - Ax2  - 


dx 


dx2 


(27) 


6  dx3  n!  dx" 


(28) 


If  we  keep  Ax  very  small  then  the  terms  smaller  than  Ax2  become  negligble.  If  equations 


(27)  and  (28)  are  added  we  can  obtain  an  expression  for  the  second  derivative: 

d2T(xo+Ax) 


T(Xo+2Ax)  +  T(Xo)  =  2T(Xo+Ax)  +  - 


dx2 


Ax2 


(29) 


d2T(Xo+Ax)  T(Xo+2Ax)  -  2T(xo+Ax)  +  T(xo) 


dx2 


Ax2 


(30) 


which  confirms  our  earlier  result. 

Thus,  as  long  as  Ax  is  kept  very  small,  the  following  are  valid  representations  and  can 


be  used  in  the  analysis  process: 

dT  T(Xo+Ax)  -  T(Xo)  AT 
dx  “ 


Ax 


Ax 


(31) 


d^T  T(xo+2Ax)  -  2T(Xo+Ax)  +  T(Xo) 
dx2 " 


(32) 
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As  stated  in  section  I,  the  configuration  to  be  analyzed  can  be  divided  into  small  but 
finite  volumes  which  are  considered  to  be  isothermal.  The  geometric  center  of  each 
subvolume  is  designated  as  a  node  and  the  node  can  be  representative  of  the  entire 
subvolume.  These  nodes  may  be  connected  to  each  other  through  thermal  resistances  and 
the  nodal  analysis  for  the  node  temperatures  lends  itself  to  solution  on  the  computer. 

The  first  law  of  thermodynamics  (the  law  of  the  conservation  of  energy)  states  that 
energy  cannot  be  created  nor  destroyed,  but  can  be  transformed  from  one  form  to  another. 
One  may  therefore  form  an  energy  balance  on  the  typical  node  shown  in  Figure  7.  For  ease 
of  explanation,  the  environment  will  be  considered  to  be  a  single  node,  node  100  (Figure 
8).  Examination  of  node  5,  shows  that  it  is  connected  to  nodes  2, 4,  6,  8,  14,  and  1(X).  If 
any  energy  is  directly  applied  to  or  removed  fiom  this  node,  we  would  also  include  it  in  the 
equation.  The  node  equation  for  node  5,  with  rate  of  heat  input  (qi),  becomes; 

q2  +  q4  +  q6  +  qs  +  qi4  +  qioo  -  qi  =  0  (33) 

where  each  one  of  the  q’s  with  a  numerical  subscript  represents  the  rate  of  heat  flow  from 
node  5  to  the  node  indicated  by  the  subscript. 

Equation  (33)  is  a  continuity  relationship  and  it  does  not  show  the  items  of  interest,  the 
node  temperatures.  But  q  =  ^  and  the  terms  in  equation  (33)  can  be  given  by: 


kA  .  k  Ax  Az  ^  , 

q2  -  r  AT  -  (T5  -  T2) 

^  Ay 

(34a) 

kA  ^  k  Ay  Az  ^  ^ 

q4-  ,  AT-  /  (T5-T4) 

(34b) 

kA  *  k  Ay  Az  ^  ^ 

q6-  j  AT-  /  (T5-T6) 

(34c) 

_  kA  ,  k  Ax  Az  ^  , 

qg-  T  AT-  (Ts-Tg) 

^  Ay 

(34d) 
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(34c) 


kA  k  Ax  Ay  ^  ^ 

qi4  =  -p  AT  = - ^  (Ts  -  Tu) 

^  Az 


kA  ^  ^  k  Ax  Ay  ^  ^  , 
qioo  =  7-  AT  = - - - (T5  -  Tioo) 

2 


(34f) 


where 

Tioo  is  on  the  top  face  of  node  5  and, 

Ax,  Ay,  Az  =  the  distance  between  nodes  in  the  x,  y,  and  z  directions  respectively 
If  we  assume  that  Ax  =  Ay  =  Az  =  1  the  node  equation  becomes: 

K(T5-T2)+K(T5-T4)+K(T5-T6)+K(T5-T8)+K(T5-Ti4)+2K(T5-Tioo)  =qi  (35) 
or 

-T2  -  T4  +  7T5  -  Tfi  -  Tg  -  Ti4  =  ^  +  2Tioo  (36) 

where  the  terms  on  the  right  hand  side  of  equations  (35)  and  (36)  should  be  known  values. 

A  node  equation  can  be  written  for  each  node  in  the  configuration  and  in  Figure  8 
which  contains  18  nodes,  18  node  equations  would  be  required  to  obtain  the  18  unknown 
node  temperatures.  To  make  maximum  use  of  the  thermal  analysis  software,  the  generation 
of  750  nodes  with  750  node  equations  similar  to  the  one  just  derived  would  be  required. 
There  are  a  number  of  methods  that  can  be  utilized  to  solve  these  equations.  It  is  not  our 
purpose  to  solve  these  equations  at  present  but  to  input  them  into  a  computer  program 
which  will  solve  them  for  us. 
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IV.  MICCKT  AND  THE  THERMAL  ANALYZER 

The  model  builder  program  which  will  become  part  of  the  thermal  analyzer  software 
package  is  MICCKT.  MICCKT  is  for  the  specific  microcircuit  configuration  shown  in 
Figure  9,  specifically,  the  microcircuit  is  covered  by  a  lid  and  supported  by  a  base.  The 
gap  between  the  chip  and  the  lid  is  assumed  to  be  air  filled.  Once  the  physical  properties  of 
the  microcircuit  are  entered  MICCKT  will  generate  an  output  file  in  the  correct  format  for 
use  in  the  thermal  analysis  software. 

A.  MICROCIRCUIT  MODELLING  PROCEDURE 

Recalling  the  idea  that  Ax,  Ay,  and  Az  must  be  kept  small  to  enable  the  use  of  the  finite 
difference  approach,  we  divide  the  microcircuit  into  subvolumes.  The  model  builder  will 
divide  the  chip  into  four  layers  with  the  top  and  second  layers  containing  four  times  as 
many  subvolumes  as  the  bottom  layers.  (Figure  9)  The  top  layer,  containing  the  circuitry, 
is  preset  at  0.{X)2  inches  (0.(X)508  cm)  and  the  thicknesses  of  each  of  the  next  three  layers 
are  one  third  the  overall  thickness  less  0.(X)2  inches.  The  reason  for  this  is  that  heat 
injection  will  occur  in  the  top  layer  and  accuracy  is  assisted  when  the  spatial  increment  is 
smaller.  In  accordance  with  the  procedures  already  discussed,  the  node  will  be  placed  in 
the  geometric  center  of  the  subvolume  and  will  be  representative  of  the  physical  properties 
of  that  subvolume. 

The  microcircuit  structure  is  made  of  pure  silicon.  This  presents  a  significant 
challenge  as  the  thermal  conductivity  of  silicon  varies  widely  with  temperature.  (Figure  10) 
In  order  to  determine  the  intemodal  conductances  it  is  necessary  to  modify  this  curve  for 
use  in  the  thermal  analyzer.  The  thermal  conductance  is  function  of  the  area  of  the  path  of 
transfer,  the  path  length,  and  the  thermal  conductivity.  Because  the  area  and  path  length  are 
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established  by  the  internodal  distances,  the  thermal  conductivity  at  a  given  temperature  is 
simply  multiplied  by  a  scalar  value  (A/L).  It  is  the  selection  of  the  thermal  ccmductivity  for 
the  desired  temperature  that  poses  a  significant  problem. 

One  of  the  features  of  the  thermal  analyzer  is  its  ability  to  generate  a  limited  facsimile 
of  a  temperature  dependent  curve  based  on  uso-  inputs.  This  temperature  dependent  curve 
may  be  specified  using  eight  straight  line  segments  between  nine  data  points.  Each  data 
point  represents  a  temperature  and  the  corresponding  conductance  for  the  physical  and 
thermal  characteristics  of  the  chip.  Because  the  physical  dimensions  of  each  subvolume 
may  change  firom  location  to  location,  several  curves  may  be  required  to  ccxi^letely  analyze 
the  configuration. 

It  was  determined  that  eight  curves  representing  the  conduction  mode  of  transfer 
within  the  structure  were  required  for  the  model.  The  top  layer  required  two  curves, 
designated  curves  1101  and  1102,  to  represent  the  north-south  and  east- west  node 
connections  and  a  third  curve,  1 103,  to  connect  the  nodes  in  layer  one  and  layer  two.  The 
second,  third,  and  fourth  layers  use  the  same  curves,  1 104  and  1 105,  for  the  north-south 
and  east-west  connections.  The  Ax  and  Ay  distances  in  the  bottom  two  layers  are  twice 
that  in  the  second  layer  but  because  either  Ax  or  Ay  will  represent  the  path  length,  the 
curves  are  identical.  Curve  1 106  is  used  to  connect  layer  two  to  layer  three  and  curve  1 107 
to  connect  layer  three  to  layer  four.  Curve  1 108  is  used  to  connect  the  fourth  layer  nodes  to 
the  heat  sink. 

The  generation  of  each  curve  requires  the  establishment  of  the  thermal  conductivity  at 
specific  temperatures.  Since  there  is  a  nine  point  limitation,  the  temperatures  were  selected 
to  correspond  to  the  expected  operating  range  of  the  chip.  The  points  selected  were: 
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Thermal  conductivity  W/cm-®C 
1.68 
1.49 


°C 

0.0 

25.0 


60.0 

1.27 

85.0 

1.14 

120.0 

1.01 

160.0 

0.88 

205.0 

0.77 

255.0 

0.66 

310.0 

0.55 

These  points  were  selected  because  they  yield  a  piecewise  continuous  curve  that  closely 
resembles  the  original  curve. 

The  thermal  interaction  between  the  top  of  the  chip  and  the  lid  of  the  container  must  be 
addressed.  The  modes  of  heat  transfer  involved  are,  conduction  and  radiation  across  the  air 
layer.  The  linearization  of  radiation  by  the  analyzer  leads  to  excessive  running  time.  Thus, 
it  would  be  beneficial  to  eliminate  the  non  linear  effect  of  the  radiation  if  at  all  possible.  A 
comparision  of  the  heat  transferred  by  conduction  to  that  transferred  by  radiation  will 
provide  additional  information  to  assist  in  justifying  an  alternate  approach  to  the  inclusion 
of  radiation. 

The  rate  of  heat  per  unit  area  for  conduction  through  the  air  gap  is: 

^lcond  =  £;(Ts-T,)  (37) 

and  for  radiation: 

^  U  =  o  Fa  Fe  (T,2  +  T,2)  (T,  +  T,)  (T,  -  T,)  (38) 

and  a  comparison  of  these  for  a  typical  case  leads  to  a  means  for  the  determination  of  the 
relative  importance  of  each  in  the  heat  transfer  process.  For  example,  using  the  values: 
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Tj  =  temperature  of  the  chip  =  400  K 
Tr  =  temperature  of  the  lid  =  300  K 
k  =  thermal  conductivity  of  air  at  350  K  »  2.97(1CH)  W/cm-K 
L  =  path  length  (air  gap  distance)  =  0.254  cm 
Fa  =  shape  arrangement  =1.0 

Fe  =  emissivity  factor,  for  parallel  plates  =  — j - ^ 

where  £c  Sl 

Ec  =  emissivity  of  the  chip  =  0.7 
El  =  emissivity  of  the  lid  =  0.8 
the  rate  of  heat  per  unit  area  for  conduction  becomes: 

^lcond=  1.17(10-1)  W/cm2 
and  for  radiation: 

^  lad  =  5.9(10-2)  W/cm2 

For  the  selected  conditions,  the  value  for  conduction  transfer  is  approximately  twice  that  for 
the  radiation  transfer.  The  difference  between  the  heat  transferred  through  the  silicon  and 
the  heat  transferred  across  the  air  layer  will  be  on  the  order  of  several  magnitudes.  In  this 
example,  for  a  chip  0.05  cm  thick  at  400  K,  the  heat  transfer  into  the  second  layer  of  the 
chip  would  be  109.75  W/cm2.  Therefore  the  vast  majority  of  heat  transfer  occurs 
downward  into  the  chip.  Due  to  it’s  poor  heat  transfer  characteristics,  air  essentially 
performs  as  an  insulator.  Because  radiation  is  even  less  effective  as  a  heat  transfer 
mechanism  the  impact  on  the  final  solution  is  minimal. 

Although  radiation  effects  are  a  small  percentage  of  the  total  heat  transferred,  these 
effects  will  not  be  neglected  in  the  model  builder.  It  will  be  incorporated  into  the  thermal 
conductivity  of  air  through  the  use  of  a  correction  factor.  The  correction  factor  will  be  a 
ratio  of  the  magnitudes  of  the  heat  per  unit  area  for  radiation  to  the  heat  per  unit  area  of 
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conduction  across  the  air  layer.  This  correction  factor  can  be  used  to  derive  a  thermal 
conductance  for  heat  flow  across  the  air  gap  at  a  given  temperature  difference  that  involves 
both  modes  of  heat  transfer 
K=^(Ratio) 

Bearing  in  mind  that  the  thermal  transfer  tu^oss  the  air  gap  by  conduction  and  radiation 
is  significantly  smaller  than  that  through  silicon,  and  the  desire  to  limit  the  running  time  of 
the  thermal  analyzer,  we  choose  to  linearize  the  radiative  transfer  mechanism  in  this 
manner.  Conduction,  with  the  correction  factor  added,  will  represent  the  entire  heat 
transfer  mechanism  between  the  chip  and  the  lid. 

Because  the  thermal  conductivity  of  air  is  temperature  dependent,  an  additional  curve 
must  be  incorporated  into  the  model  builder,  litis  curve  will  be  designated  as  curve  1 109 
and  will  be  implemented  in  the  same  fashion  as  the  silicon  conductance  curves.  The 
correction  factor  will,  of  course,  automatically  be  entered  into  this  curve. 

B.  FEATURES  OF  THE  THERMAL  ANALYZER 

The  purpose  of  the  model  builder  is  to  provide  input  data  to  the  thermal  analyzer.  This 
data  must  be  organized  in  a  format  which  is  consistent  with  the  format  of  the  thermal 
analyzer.  The  input  file  will  consist  of  five  data  lines  and  as  many  as  seven  data  sets. 
Detailed  information  regarding  the  input  data  file  may  be  obtained  in  the  user’s  manual.  In 
this  section  we  describe  features  of  specific  interest  as  they  pertain  to  the  model  builder. 
This  will  help  the  reader  gain  an  understanding  of  the  features  of  the  thermal  analyzer . 

Line  1  is  a  title  line  which  may  contain  as  many  as  79  alphanumeric  characters.  This 
title  will  be  used  to  identify  the  output  files. 

Line  2  is  the  problem  data  line  and  will  describe  inputs  pertaining  to  the  problem  at 
hand.  Of  the  items  in  this  line,  the  user  will  control  the  number  of  nodes  under 
consideration  and  the  units  used.  Other  items  identify  the  capabilities  of  the  thermal 
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analyzer  to  handle  situations  involving  heaters,  modes  of  heat  transfer  with  unique 
exponents,  nodes  which  receive  secondary  heat  inputs,  nodes  controlling  fast  heat,  and  the 
number  of  temperature  dependent  heat  input  curves.  These  entries  are  not  tq)plicable  to  this 
model  builder.  The  number  of  constant  temperature  values  and  number  of  temperature 
dependent  curves  are  applicable  and  will  be  set  by  the  model  builder. 

Line  3  is  beyond  the  user’s  control  and  is  part  of  the  thermal  analysis  software;  all 
entries  in  this  line  are  automatic. 

Line  4  is  the  problem  capability  line.  This  line  defines  the  maximiiTn  values  that  may 
be  used  for  each  item  in  line  2.  Of  particular  interest  are  the  maximum  number  of  nodes  to 
be  considered,  750,  the  maximum  number  of  temperature  dependent  curves,  15,  and  the 
maximum  number  of  constant  temperatures,  50.  The  significance  of  entering  the  750  node 
maximum  becomes  clearer  when  we  look  at  the  node  assignments  to  the  constant 
temperature  inputs.  The  first  constant  temperature  will  be  assigned  to  node  751  and 
subsequent  constant  temperatures  will  be  assigned  in  order  up  to  node  800. 

Line  5  items  place  certain  controls  on  the  computational  processes  of  the  thermal 
analyzer.  These  items  are  automatically  entered  by  the  model  builder.  Item  1  describes  the 
accuracy  between  iterations.  The  thermal  analyzer  uses  an  iterative  method  to  arrive  at  a 
solution.  As  soon  as  the  variation  of  every  temperature  computed  does  not  differ  from  its 
previously  computed  value  by  the  error  criterion  specified  here,  computation  stops.  There 
is  a  trade  off  between  accuracy  and  computer  time  involved  in  selecting  this  number.  Item 
2  is  a  damping  value  which  is  used  to  prevent  temperature  oscillations  between  iterations. 
Item  3  is  the  maximum  number  of  iterations  to  be  performed  by  the  computer.  This 
number  is  specifred  to  prevent  long  computer  runs  in  the  event  of  faulty  data.  Item  4  is  a 
convergence  factor.  Excessive  run  times  may  be  induced  by  the  damping  factors  if  there  is 
a  slow  convergence  of  either  overdamped  or  underdamped  temperature  values.  If  it  has 
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been  determined  that  this  is  occurring,  the  convergence  factor  will  increase  the  convergence 
due  to  damping  in  order  to  reduce  computer  time.  Item  S  specifies  the  initial  temperature 
estimate  fcH*  use  in  the  iteration  process. 

Input  data  set  1  contains  the  temperature  dependent  coefficients  for  as  many  as  IS 
curves.  These  curves  and  their  structure  have  been  previously  described  in  the  microcircuit 
modelling  sectitm. 

Input  data  set  2  contains  up  to  50  constant  temperatures. 

Input  data  set  3  involves  heaters  for  fast  warm  up  and  is  not  used  in  this  model. 

Input  data  set  4  contains  the  information  for  the  construction  of  the  n  node  equations. 
For  each  node  there  will  be  two  lines  of  data.  The  odd  numbered  lines  specify  which 
nodes  are  connected  to  the  node  in  question  and  the  mode  of  heat  transfer.  Each  mode  of 
heat  transfer  is  assigned  a  numerical  value  from  1  to  10.  Our  model  builder  is  concerned 
only  with  conduction  which  has  been  designated  by  the  number  one  (1).  For  example,  if 
node  25  was  connected  to  node  2  by  conduction,  in  line  3  of  data  set  4  an  entry  of  251 
would  appear.  The  even  numbered  lines  specify  the  intemodal  conductances.  For  the 
model  builder  this  line  will  contain  the  curve  number  specified  in  data  set  1  corresponding 
to  the  connecting  nodes. 

Input  data  set  5  specifies  the  use  of  unique  exponents  and  is  not  used  in  this  model. 

Input  data  set  6  specifies  the  initial  temperature  guesses  corresponding  to  the  number 
of  nodes  receiving  secondary  heat  and  is  not  used  in  this  model.  All  initial  temperature 
guesses  are  equal  and  are  set  in  line  5. 

Input  data  set  7  specifies  the  number  of  temperature  dependent  heat  input  curves  and  is 
not  used  in  this  model. 
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C.  FEATURES  OF  MICCKT 

The  MICCKT  program  will  generate  a  file  which  can  be  read  by  the  thermal  analysis 
software  in  order  to  provide  the  thermal  analysis.  The  program  is  user  friendly  and  will 
save  numerous  man  hours  in  the  analysis  process.  Through  the  use  of  MICCKT  reaching 
the  750  node  capability  of  the  analyzer  is  now  possible  with  a  minimum  of  effort. 

The  program  is  interactive  and  will  provide  prompts  for  the  required  inputs.  After  the 
prompt  for  a  title  for  the  output  file  there  will  be  a  list  of  options  to  produce  the  number  of 
nodes  for  analysis.  These  options  are: 

1)  16  by  16,  (640  nodes) 

2)  12  by  24,  (720  nodes) 

3)  10  by  30,  (750  nodes) 

4)  8  by  32,  (640  nodes) 

5)  12  by  18,  (540  nodes) 

The  options  were  selected  to  provide  the  user  with  several  aspect  ratios  to  approximate  the 
microcircuit  structure.  Option  1  provides  the  user  with  a  1:1  aspect  ratio  while  options  2 
through  5  provide  the  aspect  ratios  of  1:2, 1:3, 1:4,  and  2:3  respectively. 

After  selecting  the  number  of  nodes  for  analysis  the  user  must  specify  the  use  of  either 
English  or  SI  units.  Once  this  selection  is  made  all  subsequent  inputs  will  conform  to  the 
chosen  units.  The  next  several  inputs  are  concerned  with  the  physical  size  of  the  chip  and 
with  initial  conditions.  The  user  will  input  the  length,  width,  and  overall  thickness  of  the 
chip,  the  initial  temperatures  of  the  lid  and  sink,  and  the  air  gap  distance  between  the  chip 
and  the  lid.  As  a  reminder,  the  user  is  limited  to  an  air  gap  thickness  of  0.254  centimeters 
(0.10  inches)  maximum. 

It  is  at  this  point  that  any  heat  inputs  to  the  chip  may  be  entered.  The  user  may  select 
from  four  options.  The  first  option  is  total  rate  of  heat  applied  to  the  chip.  The  total  rate  of 
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heat  will  be  divided  by  the  number  of  nodes  in  the  top  layer  and  evenly  distributed  to  each 
node  as  an  average  input.  The  second  option  is  average  rate  of  heat  per  unit  area,  this  will 
be  either  in  Btu/hr-in2  or  Watts/cm^.  This  average  rate  of  heat  per  area  will  be  multiplied 
the  area  of  the  chip  and  evenly  distributed  to  each  node.  The  third  option  allows  the  user  to 
select  specific  nodes  for  individual  heat  rate  inputs.  The  fourth  option  is  for  no  external 
heat  input  to  the  chip. 

Once  these  selections  are  made  a  summary  of  the  user’s  inputs  is  displayed  for 
verification.  The  user  can  review  these  inputs  and  make  necessary  changes  before 
continuing.  If  all  entries  are  satisfactory  the  run  will  continue  and  produce  an  ou^ut  file. 

This  output  file  will  be  in  the  the  proper  format  for  use  in  the  thermal  analyzer.  All 
data  input  lines  and  data  input  sets  will  be  automatically  generated.  In  particular,  for  input 
data  set  4,  each  node  will  be  connected  to  its  neighboring  nodes  by  the  mode  of  heat 
transfer  and  appropiate  temperature  dependent  curve. 
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C  APPENDIX  A 

C 

C 

C  THIS  PROGRAM  15  DESIGNED  TO  PROVIDE  AN  OUTPUT  FILE  IIHICH  MILL  BE 
C  READ  BY  THE  THERMAL  ANALYSIS  SOFTIIARE.  IT  IS  A  MODEL  BUILDER  FOR  A 
C  MICROCIRCUIT  STRUCTURE  OF  PURE  SILICON  AND  HILL  TAKE  INTO 
C  CONSIDERATION  IMF.  VARIATIONS  IN  THERMAL  CONDUCTIVITY  EXHIBITED  BY 
C  SILICON.  THE  STRUCTURE  IS  DIVIDED  INTO  FOUR  LAYERS  IIITH  THE  TOP  TITO 
C  LAYERS  CONTAINING  FOUR  TIMES  AS  MANY  NODES  AS  THE  BOTTOM  TMO  LAYERS, 

C  THE  STRUCTURE  IS  COVERED  IIITH  A  LID  AND  THE  GAP  BETMEEN  THE  LID  AND 
C  CHIP  IS  AIR  FILLED. 

C  COMMENTS  ARE  ENTERED  THROUGHOUT  THIS  PROGRAM  IN  ORDER  TO  EXPLAIN 
C  THE  DEFINITION  OF  SYMBOLS  USED  OR  TO  EXPLAIN  THE  LOGIC  INVOLVED. 

C  THIS  PROGRAM  MAS  DEVELOPED  AS  PART  OF  A  THESIS  BY  LCDR  JOSEPH  A. 

C  MILLHELM,  "COMPUTER  AIDED  THERMAL  ANALYSIS  OF  A  MICROCIRCUIT 
C  STRUCTURE",  DEC.  1990.  FURTHER  INFORMATION  CAN  BE  OBTAINED  IN  THIS 
C  THESIS. 

C 

C 

CHARACTER  *  79  TITLE 
CHARACTER  x  1  RESP,  YES,  ANS,  CHANS 
C 

C  THE  TITLE  ALLOWS  FOR  UP  TO  79  ALPHANUMERIC  CHARACTERS  TO  BE 
C  ENTERED  TO  DESCRIBE  THE  OUTPUT  FILE.  RESP  IS  THE  RESPONSE  TO  MAKING 
C  ANY  CHANGES  TO  THE  SUMMARY  WHEN  IT  APPEARS.  YES  IS  THE  RESPOtlSE  TO 
C  DIFFERENT  QUESTIONS  AS  THEY  OCCUR,  ANY  OTHER  RESPONSE  WILL  RESULT  IN 
C  CONTINUATION.  ANS  IS  THE  RESPONSE  TO  MAKING  ANY  CHANGES  TO  THE  HEAT 
C  INPUTS.  CHANS  IS  THE  RESPONSE  TO  ALL  CHANGES  MADE  IN  THE  HEAT  INPUTS. 
C 
C 

INTEGER  R,  COUNT,  CNTR,  NODE,  I.  J.  XROW.  YROW,  SELECT,  NDE, 
+CHNODE,  LVL,  COUNTER,  NX,  NY,  NR,  NL ,  CTR,  KNTR, 

+ALPHA,  BETA,  G,^MMA,  DELTA.  NP.  N 
C 

INTr(N  R  CONST.  MIR,  rxP.  SFCHT  ,  COn  CV,  FASTIU,  HICV, 

4C0NIK()L1.  CNMIRULiJ,  CONIRlILi,  NODLS.  KMAX.  HMAX,  XP,  SHT,  EFCV, 
+TCV,  STHT,  MAXIIS,  DMM,  UNK,  UNIT,  LIDN,  SNKN 


THE  FOLLOWING  ARE  A  LIST  OF  THE  ABBREVIATIONS  USED  IN  THE  INTEGER 
LISTINGi 

COUNT,  CNTR,  COUNTER,  CTR,  AND  KNTR  =  COUNTERS 
NODE,  NDE  =  NODE  NUMBERS  OF  THE  STRUCTURE 
I,  J,  AND  N  =  USED  IN  DO  LOOPS,  COUNTERS 
XROW  =  NUMBER  OF  NODES  IN  THE  X  DIRECTION 
YROW  =  NUMBER  OF  NODES  IN  THE  Y  DIRECTION 
NX  =  XROW/2.  NUMBER  OF  NODES  IN  THE  X  DIRECTION  OF  THE 
BOTTOM  TWO  LAYERS 

NY  =  YROW/2,  NUMBER  OF  NODES  IN  THE  Y  DIRECTION  OF  THE 
BOTTOM  TWO  LAYERS 

NR  =  NEXT  ROW,  TO  DETERMINE  THE  NODE  ADJACENT  TO  THE 
NODE  UNULR  CONSIDERATION,  +  OR  -  Y  DIRECTION 
NL  =  NEXT  LAYER,  TO  DETERMINE  THE  NODE  ADJACENT  TO  THE 
NODE  UNDER  CONSIDERATION,  ABOVE  OR  BELOW 
R  =  TOTAL  NUMBER  OF  NODES  SELECTED 
NP  =  NEXT  POINT,  COUNTER 

X,  Y,  AND  Z  =  THREE  DIMENSION  COORDINATE  SYSTEM 
SELECT  =  NUMBER  OF  THE  OPTION  SELECTED  FOR  NUMBER  OF 
NODES  TO  BE  USED 

CHNODE  =  CHANGE  NODE.  TO  CHANGE  THE  NODE  SELECTED  FOR 
HEAT  INPUT 

LVL  =  LEVEL,  If*  THIS  CASE  SET  TO  A,  ALTHOUGH  HOT 
SPECIFICALLY  USED  IN  THIS  PROGRAM  IF  FUTURE 
APPLICATIONS  REQUIRE  A  VARIABLE  LAYER  SELECTION  THE 
BASIS  FOR  THIS  PURPOSE  ARE  IMBEDDED  IN  THIS  PROGRAM 
ALPHA,  BETA,  GAMMA,  AND  DELTA  =  ALGORITHMS  USED  TO 

COtINECT  THE  NODES  IN  LAYER  TMO  TO  THE  tfODES  Iff  LAYER 
THREE  due  to  THE  A.l  RATIO  BEIMLLN  THE  LAYERS 
UNIT  =  RESPONSE  TO  51  OR  ENGLISH  UNITS 
SNKN  =  NUMBER  OF  SINK  NODES 

LIDN  =  NUMBER  OF  NODES  TO  THE  LID  AND  AIR  GAP 
ALL  OTHER  INTEGER  ENTRIES  INVOLVE  NECESSARY  INPUTS  FOR  THE 
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C  THERMAL  AtlALYZER  AtID  ARE  USED  HERE  FOR  ACCOUHTItIG  PURPOSES. 

C 

C 

INTEGER  8(750,800  ).  TMP(IO),  KOUNTER( 7 50  ) ,  SUB1(9) 

C 

C 

C  B  ARRAY  =  »IOnE  NUMBER  ARRAY 

C  TMP  AND  SUBl  ARRAYS  =  TEMPORARY  ARRAYS  USED  FOR  COUNTING 

C  KOUNTER  ARRAY  =  USED  AS  A  COUNTER 

C 
C 

REAL  MATTS.  TOTHT,  HTPAR.  AVG.  AVGl ,  lENGTH,  WIDTH,  EP5L ,  EPSD, 
+EMS,  TX,  TAV.  HGAP.  HRAD,  FE,  LID.  SINK,  CHOICE,  LINE,  CHWATT, 

+  THICK,  THl,  TH2.  DEL'X,  OELY,  DC121,  DEIZ2,  RAT,  FACTOR, 

+DELZ3,  DELZA,  GAP.  TEMPI,  TEMP2,  TEMPS,  TEHP4,  TEMPS,  TEMP6, 
+TEMP7,  TEMPS,  TEMP9,  ANSW,  KAIR,  LIDK, 

+K1,  K2,  K3,  KA,  K5.  K6 ,  K7,  K8,  ACC,  DAMP,  CONVFAC,  INTEMP 


MATTS  =  HEAT  INPUT  TO  THE  MODES,  CAN  REPRESENT  HEAT  INPUT 
IN  EITHER  MATTS  OR  BTU/HR 
TOTHT  =  TOTAL  HEAT  APPLIED  TO  THE  CHIP 
HTPAR  =  HEAT  PER  UNIT  AREA  APPLIED  TO  THE  CHIP 
AVG  =  TOTAL  HEAT  DIVIDED  BY  THE  NUMBER  OF  NODES  IN  THE 
FIRST  layer,  PROVIDES  EQUAL  DISTRIBUTION  OF  HEAT 
AVGl  =  HEAT  PER  AREA  TIMES  THE  AREA  DIVIDED  BY  THE  NUMBER 
OF  NODES  IN  THE  FIRST  LAYER 
LENGTH  =  lENGTH  OF  THE  CHIP,  X  DIRECTION 
MIDIH  =  MIDIH  OF  THE  CHIP,  Y  DIRECTION 
EPSL  =  EMISSIVITY  OF  THE  LID 
EPSD  =  EMISSIVITY  OF  THE  CHIP  (DIE) 

EMS  =  EMISSIVITY  EQUIVALENCE,  1/EPSD  +  1/EPSL  -  1 
FE  =  INVERSE  OF  EMS,  EMISSIVITY  FACIOR 
HRAI)  -  KADI  AI  ION  CONSIANT  OF  PROPOR 11  ONAl  I  T  Y  FOR  AIR 
HGAP  =  COMDIK  I  ION  CONSIANT  OF  PROPOR I  1 ONAL I  I Y  FOR  AIR 
RAI  =  RAIIC)  OF  HRAD  TO  HGAP,  CORRECTION  FACTOR 
FACTOR  =  HRAT,  MULTIPLIER  TO  PUT  CORRECTION  INTO  THE 
THERMAL  CONDUCTIVITY  CURVE  FOR  AIR 
TX  =  STARIING  TEMPERATURE  TO  CALCULATE  THERMAL 
CONDUCIIVITY  OF  AIR  TO  GENERATE  CURVE  1109 
TAV  =  AVERAGE  TEMPERATURE  OF  THE  LID  AND  TX 
K1  THROUGH  K8  =  THERMAL  CONDUCTIVITY  CURVES  OF  SILICON 
REFER  TO  THESIS  FOR  MORE  DETAILED  EXPLANATION 
KAIR  =  THERMAL  CONDUCTIVITY  CURVE  OF  AIR  (1109) 

LID  =  TEMPERATURE  OF  THE  LID 
SINK  =  TEMPERATURE  OF  THE  SINK 
CHOICE  =  CHOICE  OF  THE  METHOD  OF  HEAT  INPUT 
LINE  =  LINE  NUMBER  THAT  MILL  BE  CHANGED  AFTER  SUMMARY 
CHMATT  =  CHANGE  THE  VALUE  OF  THE  HEAT  INPUT 
THICK  =  OVERALL  THICKNESS  OF  THE  CHIP 
THl  =  THICKNESS  OF  THE  TOP  LAYER,  PRESET  AT  0.002  IN 
(0.00508  CM) 

TH2  =  THICKNESS  OF  THE  SECOND,  THIRD,  AND  FOURTH  LAYERS 
(THICK-THl  )/3 

DELX  =  LENGTH  DIVIDED  BY  THE  NUMBER  OF  NODES  IN  THE  X 
DIRECTION  TO  DETERMINE  SEPARATION  DISTANCE 
DELY  =  MIDIH  DIVIDED  BY  THE  NUMBER  OF  NODES  IN  THE  Y 
DIRECTION  TO  DETERMINE  SEPARATION  DISTANCE 
DELZl  =  NODE  SEPARATION  IN  THE  Z  DIRECTION,  LEVEL  1 

DELZ2  =  NODE  SEPARATION  IN  THE  Z  DIRECTION,  LEVELS  2-A 

DELZ3  =  NODE  SEPARATION  IN  THE  Z  DIRECTION,  BEliTCEN 

LAYERS  2  AND  3 

DELZA  =  NODE  SEPARATION  IN  THE  Z  DIRECTION,  BET5IEEN 
LAYER  A  AND  THE  SINK 

GAP  =  DISTANCE  BETWEEN  THE  TOP  OF  THE  CHIP  AND  THE  LID 
PRESET  TO  A  MAXIMUM  OF  0.10  IN  {0.25A  CM) 

TEMPI  THROUGH  TEMP9  =  TEMPERATURES  SEIFCTEI)  FROM  THE 

SILICON  THERMAL  COtMlUC T I V I T Y  CURVE  TO  ESTABLISH  THE 
CURVES  FOR  USE  IN  THE  THERMAL  ANALYZER 
ANSM  =  RESPONSE  TO  THE  OPTIONS  OFFERED  FOR  NUMBER  OF 
NODES  TO  BE  USED 

LIDK  =  LID  TEMPERATURE  IN  KELVIN 
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AIL  OTHFR  RIAL  EtITRIES  IIIVOIVE  UECESr.ARY  IMPIITS  FOR  IMF 
IIII.KMAL  ANAIYZEK  AMD  ARE  USED  HLKL  FOR  ACCmiMIlMO  I’MRI’IISES 


REAL  A(750,8n0),  TEMP(IO).  L10T(50),  SiriKT(<*), 
4THC0MU(9),  CURVE(8,9),  C(8),  SUB(9),  TMA1R(9),  KAR(9) 


A  ARRAY  =  HEAT  INPUT  ARRAY 

TEMP  AND  SUB  ARRAYS  =  TEMPORARY  ARRAYS  FOR  COUNTING 
LIDT  ARRAY  =  INPUTS  CONSTANT  TEMPERATURE  OF  THE  LID 
FOR  THE  THERMAL  ANALYZER  BASED  ON  REQUIREMENTS 
SINKT  ARRAY  =  INPUTS  CONSTANT  TEMPERATURE  OF  THE  SINK 
FOR  THE  THERMAL  ANALYZER  BASED  ON  REQUIREMENTS 
THCOND  ARRAY  =  THERMAL  CONDUCIVITY  OF  SILICON  AT  THE 
TEMPERATURES  SELECTED  IN  TEMPI  -  TEMP9 
C  ARRAY  =  AREA  DIVIDED  BY  THE  PATH  LENGTH  BASED  ON  INPUTS 
CURVE  ARRAY  =  C  ARRAY  TIMES  THE  THCOND  ARRARY 

THE  TEMPI  -  TEMP9  POINTS  IN  CONJUNCTION  WITH  THE 
CURVE  ARRAY  PROVIDE  THE  THERMAL  ANALYZER  WITH  THE 
DATA  TO  GENERATE  THE  1101  -  1108  CURVES 
•THAIR  ARRAY  =  EQUATION  TO  GENERATE  THERMAL  CONDUCTIVITY 
OF  AIR  GIVEN  A  TEMPERATURE 
KAR  ARRAY  =  THE  K  VALUES  FOR  THE  CURVE  1109 


PARAMETER  (YES= 'Y' ) 


PRINT  K, 'ENTER  THE  NAME  OF  YOUR  FILE.' 

READ  !«.  TITLE 

PRINT  X.'YOU  MAY  SELECT  ONE  OF  THE  FOLLOWING  OPTIONS.' 

30  PR  I II I  16  BY  16,  Y I  FIDS  660  NODFS' 

PRINI  K.’Z.)  IZ  Ur  26.  YIELDS  720  NODES' 

PRIlir  K.'3.)  10  BY  30,  YIELDS  750  NODES' 

PRIlir  K.'6.)  8  BY  32.  YIELDS  660  NODES' 

PRItir  X,'5.)  12  BY  18,  YIELDS  560  NODES' 

READ  X,  AllSU 
C 
C 

C  THESE  OPTIONS  WILL  ALLOW  THE  USER  TO  SELECT  THE  NUMBER  OF  NODES 
C  FOR  ANALYSIS.  THESE  OPTIONS  WERE  DESIGNED  TO  ALLOW  THE  USER  TO 
C  SELECT  FROM  A  LENGTH  TO  WIDTH  RATIO  OF  1:1,  1:2,  1:3,  1:6,  OR  2:3. 
C 
C 

IF  ( AN5W.NE. 1 . AND.AN5W.NE.2,AND.ANSW.NE.3.AND. ANSW.NE.6 . AND. 
+ANSW.NE.5)  THEN 

PRINT  x.'ONLY  OPTIONS  1 . , 2 . , 3 . , 6 . , OR  5.  ARE  ALLOWED.' 

GOTO  30 

ELSCIF  (ANSW.EQ.l)  THEN 
XROW=16 
YROW=16 
C0NST=18 

ELSEIF  (ANSW.EQ.Z)  THEN 
XR(ni=12 
YR()W  =  26 
CON5T=20 

ELSEIF  (ANSW.Eq.3)  THEN 
XROW=10 
YRnW=30 
COH5T=20 

ELSCIF  (AN3W.EQ.6)  THEN 
XROW=8 
YROW=32 
CnN5T=18 

ELSI  II  (Ali;;il.  FQ.  5)  THEN 
XR()ll=li! 

YR()II=18 
■COIIST  =  15 
ENDIF 


36 


ooooooooo 


10  PRINT  X, 'ENTER  1.)  FOR  SI  OR  2.)  FOR  ENGLISH 
+UNI15. ' 

READ  X,  UNIT 


THE  UNITS  USED  FOR  INPUTS  IN  THIS  PROGRAM  Hill  DEt 
SI;  CENI  IMF.HRS,  HATTS,  UEO  CEIIIIGRAUE 
ENGLISH;  INCHES,  BTU/HR,  DEG  FAHRENHEIT 
OTHER  PARAMETERS  USED  IN  THIS  PROGRAM  HILL  AUTOMATICALLY  CONFORM 
TO  THE  CHOSEN  UNITS. 


C 


C 

C 


C 


C 

120 


C 


C 


C 


C 


C 

C 

c 


c 


IF  (UNIT.NE. 1 .AND.UNIT.NE.2)THEN 
PRINT  X, 'TRY  AGAIN* 

GOTO  10 
ELSE 

CONTINUE 

EHDIF 


IF  (UNIT.EQ.l)  THEN 

PRINT  X,' INPUT  CHIP  DIMENSIONS  IN  CENTIMETERS' 
ELSE 

PRINT  X,' INPUT  CHIP  DIMENSIONS  IN  INCHES’ 

ENDIF 


PRINT  X, 'ENTER  LENGTHi  ' 
READ  X,  LENGTH 


PRINT  X, 'ENTER  MIDTH:  ' 

READ  X,  MIDTH 

PRINT  X,'EflTER  THICKNESS:’ 

READ  X,  THICK 

PRINT  X, 'ENTER  THE  THICKNESS  OF  THE  AIR  LAYER  DETMEEN  THE  LID’ 
PRINT  X.'AND  THE  TOP  OF  THE  DIE' ’ 

READ  X,  GAP 


IF  (UNIT.EQ.l .AND. GAP. OT.O.eSA)  THEN 

PRINT  x,'AIR  LAYER  THICKNESS  EXCEEDS  MAXIMUM  ALLOMED  (0.25ACM)’ 
GO  IN  120 

ELSEIF  (UNIT. EQ. 2. AND. GAP. GT. 0.10)  THEN 

PRINT  x,'AIR  LAYER  THICKNESS  EXCEEDS  MAXIMUM  ALLOMED  (0.10  IN)’ 
GOTO  120 
ELSE 

CONTINUE 

ENDIF 

PRINT  X, 'ENTER  THE  TEMPERATURE  OF  THE  LID.’ 

READ  X,  LID 

PRINT  X, 'ENTER  THE  TEMPERATURE  OF  THE  SINK.’ 

READ  X,  SINK 

PRINT  X, 'ENTER  THE  EMISSIVITY  OF  THE  LID.’ 

READ  X,  EPSL 

PRINT  X, 'ENTER  THE  EMISSIVITY  OF  THE  DIE.’ 

READ  X,  EPSD 


NX  =  XRnil/2 

NY  =  YR()M/2 

LVL  =  'i 

NR=XROM 

NL=XROHXYRnH 

R  =  xr<OII><YROMX2 . 5 

COUNl LR=1 

DO  26  1  =  1, R 
DO  27  J=l,800 
A(  I , J  )  =0 . 0 
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n( I , j)=o 
('.util  I  null 
coin  IIIUI: 


ItllTlALIZCS  THH  MATRICES  WITH  ZEROES.  THIS  15  DONE  AS  THE  SIZE  OF 
THESE  MAIRICES  COULD  CAUSE  EXCESS  COMPUTER  TIME  IP  HOI  DONE  111  THIS 
FASHIOII. 


0  DO  11  1  =  1,  R 

A(I,770)=0.0 
n( 1,7701=0 
COHTIHUE 


III  THE  EVENT  THE  USER  WISHES  TO  MAKE  SIGHIFICAHT  CHANGES  TO  THE 
NODE  HEAT  INPUTS  THIS  DO  LOOP  HILL  RESET  THE  MODES  AND  HEAT  INPUTS 
TO  ZERO.  THIS  IS  DONE  IN  THE  EVENT  THAT  THE  MODE  OF  APPLICATION  HAS 
CHANCED  UR  IF  ONE  OR  MURE  NODES  ARE  TO  BE  ELIMINATED. 


IF  (UNIT.EQ.l)  THEN 

PRINT  x.'HEAT  INPUTS  WILL  BE  IN  WATTS.' 

TH1=0. 00503 
INTEnP=100 
ELSE 

PRINT  x.'HEAT  INPUTS  WILL  BE  IN  BTU/HR.' 

Till  =0 . 002 
ItlTEMP  =  212 
ENDIF 
C 
C 

C  BASED  ON  THE  UNITS  SELECTED  THIS  STATEMENT  SERVES  AS  A  PROMPT  AND 

C  TO  SET  HIE  THICKNESS  01  THE  TOP  LAYER  AND  INITIAL  TEMPERATURE  TOR 
C  THE  THERMAL  ANALYZER  IN  THE  PROPER  UNITS. 

C 

C 

PRINT  K,'HEAT  INPUT  TO  THE  MICROCIRCUIT  CAN  BE  ENTERED  IN  THE' 
PRIIIl  X,  M  01.1  LOIHIKT  MANNER.  CHOOSE  1.)  TO  ENILR  AS  TOIAL  HEAT' 
PRItll  X, 'AI'I’LIED  TO  THE  CHIP,  2.)  TO  ENTER  AS  THE  AVERAGE  HEAT' 
PRINT  x.'PER  AREA,  3.)  TO  ENTER  NODE  BY  NODE,  OR  A.)  FOR  NO  HEAT', 
PRINT  X,  'INPUT. ' 

READ  X,  CHOICE 
C 

IF  (CHOICE. EQ.n  THEN 

PRINT  X, 'ENTER  TOTAL  HEAT  APPLIED  TO  THE  CHIP' 

RFiAD  X,  TOTHT 
AVG=TOTHT/NL 
DO  12  N=1,NL 
A(N,770)=AVG 
0(11,7701=9991 

12  CONTINUE 

ELSEIF  (CHOICE. EQ. 21  THEN 

PRINT  X, 'ENTER  AVERAGE  HEAT  PER  UNIT  AREA' 

READ  X,  HI  PAR 

AVGI=(HTPARxLCNGTHXUIDTH1/NL 
DO  13  N=1,NL 

A(N,7701=AVG1 
D(N, 7701=9991 

13  CONTINUE 
C 

C 

C  THE  NUMBER  9991  IS  RECOGNIZED  BY  THE  THERMAL  ANALYZER  TO  BE  A 
C  CONDUCTIVE  HEAT  INPUT  FROM  AN  EXTERNAL  SOURCE.  THE  NODE  NUMBER  770 
C  IS  ASSIGNED  TO  TRACK  THIS  HEAT  INPUT. 

C 

C 

ELSEIF  (CHOICE. EQ. 31  THEN 

PRINT  x,'YOU  MAY  SELECT  HEAT  INPUTS  TO  THE  INDIVIDUAL  NODES.' 

Ill  PRINT  x,'YOU  ARE  LIMITED  TO  A  MAXIMUM  OF  256  INPUTS  FOR  THE' 
PRINT  x,'6A0  NUDE  MODELS,  288  INPUTS  FOR  THE  720  NODE  MODEL,' 
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PRINT  K,'300  INPUTS  FOR  THE  750  NODE  MODEL,  OR  216  INPUTS' 
PRINI  K.'FOR  1ME  5<^i0  NODE  MODEL.' 

PRIIM  K.'SElECr  INE  NUMBER  OF  NODES  THAT  HAVE  HEAT  INPUTS.' 
READ  X,  SELECT 


THE  REASON  FOR  THESE  OPTIONS  IS  THAT  THEY  REPRESENT  THE  NUtinER  OF 
NODES  IN  THE  TOP  LAYER.  THESE  ARE  THE  ONLY  NODES  ALLOHED  TO  HAVE  AN 
EXTERNAL  HEAT  INPUT  BY  THIS  PROGRAM. 


IF  (SELECT. GT.NL)  THEN 

PRINT  X, 'NUMBER  OF  NODES  SELECTED  EXCEEDS  MAXIMUM  ALLOHED.' 
GOTO  111 
ELSE 

CONTINUE 

ENDIF 

DO  lA  1=1, SELECT 

PRINT  X, 'SELECT  THE  NODE  TO  BE  INPUTTED.' 

READ  X,  NODE 

PRINT  X, 'ENTER  THE  HEAT  INPUT  AT  THAT  NODE.' 

READ  X,  HATTS 
A(N0nE,770)=HATTS 
B(NODE, 7701=9991 
KOUNTER( I )=I 


lA 

C 


17 

C 

112 


C 


CONTINUE 

ELSEIC  (CHOICE. EQ. A)  THEN 
GOTO  17 

ELSEIF  ( CHO ICE. NE. 1. AND. CHOICE. NE. 2. AND. CHOICE. NE. 3. AND. 
4CH0ICE.NE.A)  THEN 

PRINT  X, 'TRY  AGAIN' 

GO  10  110 
ENDir 
CONTINUE 


PRINT  X, 'SUMMARY  OF 
IF  (ONIT.EO.l)  IHLH 
PRINI  '(A.GIZ.A) ' , • 
+LENGIH 

PRINT  '(A,G12.A)',  • 
+HIDTH 

PRINT  '(A,G12.A)',  ' 
+THICK 


INPUT  DATA.' 

1. )  LENGTH  OF  THE  CHIP  IS  (CM)i 

2. )  HIDTH  OF  THE  CHIP  IS  (CM): 

3. )  THICKNESS  OF  THE  CHIP  IS  (CM): 


PRINT  '(A,G12.A)', '  A.)  THICKNESS  OF  THE  AIR  LAYER  (CM)i 
+GAP 

PRINT  '(A,G12.A)','  5.)  TEMPERATURE  OF  THE  LID  IS  (DEG  C): 

+  LID 

PRINT  '(A,G12.A)', '  6.)  TEMPERATURE  OF  THE  SINK  IS  (DEG  C)i 
+SINK 

PRINT  '(A,G12.A)', •  7.)  THE  EMISSIVITY  OF  THE  LID  IS: 

+  EPSL 

PRINT  '(A,G12.A)», •  8.)  THE  EMISSIVITY  OF  THE  DIE  IS: 

+  EPSD 

PRINT  x,'9.)  THE  HEAT  INPUTS  TO  THE  NODES  ARE:' 


IF  (CHOICE. EQ. 1 )  THEN 

PRINT  '(A,G12.A,A)', •  AVERAGE  HEAT  APPLIED  TO  EACH  NODE  IS:  », 
4AVG, •  MATTS' 

ELSEIF  (CHOICE. EQ. 2)  THEN 

PRINI  • (A,U12.A,A) • , '  HEAT  APPLIED  TO  EACH  NODE  ISi  *, 

4AVG1,'  MATTS' 

ELSEIF  (CHOICE. EQ. 3)  THEN 
DO  15  1=1, NL 

IF  (B(I,770).EQ.0)  THEN 
GOTO  15 
ELSE 
NDE  =  I 

PRINT  '(A, 13)',  •  NODE  '.NDE 
PRINT  '(A,G12.A,/)',  '  MATTS ' , A( I , 77 0 ) 
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cot)  I  nine 

EL  sen  (ciioice.rQ.^i)  then 

PRINT  i(,'|IO  HEAT  INPUTS  HERE  DESIGNATED.' 

ENDIF 
ELSE 

PRINT  •(A,G12.A)', •  1.)  LENGTH  OF  THE  CHIP  IS  (IN)i 
+LENGIH 

PRINT  •(A,G12,‘i)',  '  2.)  HIDTH  OF  THE  CHIP  IS  (IN)i 
+UIDTH 

PRINT  •(A.GI2.<*)'.  •  3.)  THICKNESS  OF  THE  CHIP  IS  {IN)i 
+THICK 

PRINT  •(A,G12.‘i)',  '  <« .  )  THICKNESS  OF  THE  AIR  LAYER  (IN)i 
+GAP 

PRINT  '(A.GU.A)',  •  5.)  TEMPERATURE  OF  THE  LID  IS  (DEG  F)  i 
4LID 

PRINT  '(A,G12.<i)',  '  6.)  TEMPERATURE  OF  THE  SINK  IS  (DEG  F)t 
+SINK 

PRINT  •(A,G12.'i)',  •  7.)  THE  EMISSIVITY  OF  THE  LID  ISi 
4EPSL 

PRINT  '(A,G12.A)', •  8.)  THE  EMISSIVITY  OF  THE  DIE  ISi 
4EPSD 

PRINT  *,'9.)  THE  HEAT  INPUTS  TO  THE  NODES  AREi ' 

IF  (CHOICE. EQ.l)  THEN 

PRINT  '(A.GU.A.A)','  AVERAGE  HEAT  APPLIED  TO  EACH  NODE  IS« 
4AVG, '  DTU/HR' 

ELSEIF  (CHOICE. EQ. 2)  THEN 

PRINT  •(A.G12.<i,A)',  •  HEAT  APPLIED  TO  EACH  NODE  ISt 
4AVG1.'  BTU/HR' 

ELSEIF  (CHOICE. EQ. 3)  THEN 
DO  22  1=1, NL 

IF  (B(I,770).EQ.0)  THEN 
GOlO  22 
ELSE 
NUE=I 

PRINT  '(A, 13)',  •  NODE  '/NDE 
PRINT  '(A,G12.A,/)',  '  BTU/HR • , A(  1 , 77 0 ) 

END  IF 
22  CONTI  NIIE 

ELSEIF  (CHOICE.  EO.'i)  THEN 

PRINT  x,'NO  HEAT  INPUTS  HERE  DESIGNATED.’ 

ENDIF 

C 

ENDIF 

C 

PRINT  K,'DO  YOU  NISH  TO  HAKE  ANY  CHANGES,  Y  OR  N? ' 

READ  (x.S)  RE5P 
FORMAT  (AI) 

IF  (RESP.EQ.YES)  THEN 

PRINT  K,M1UICH  LINE  ITEM  DO  YOU  WISH  TO  CHANGE?' 

READ  X,  LINE 
ELSE 

GOTO  113 
ENDIF 

IF  (LINE. EQ.l)  THEN 

PRINT  X, 'ENTER  THE  LENGTH  OF  THE  CHIP.' 

READ  X, LENGTH 
GOTO  112 

ELSEIF  (LINE.EQ.2)  THEM 

PRINT  X, 'ENTER  THE  HIDTH  OF  THE  CHIP.’ 

READ  X, HIDTH 
GOTO  112 

ELSEIF  (LINE.Eq.3)  THEN 

PRINT  X, 'ENTER  THE  THICKNESS  OF  THE  CHIP.' 

READ  X, THICK 
GOTO  112 
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c 


c 


C 


C 


C 

4 

114 


6 


C 


fiir.rip  (iniP.rQ.'i)  them 

I'KIMI  )i,<IMtEK  DIE  THICKNESS  OP  THE  AIR  LAYER.* 

READ  X,GAI> 

GOTO  112 

ELSCIF  (LIME.EQ.5)  THEM 

TRIMT  K, 'ENTER  THE  LID  TEMPERATURE.’ 

HEAD  K.LIU 
GOTO  112 

ELSEIF  (LIME.EQ.6)  THEM 

PRINT  X. 'ENTER  THE  TEMPERATURE  OF  THE  SINK.' 

READ  x.SIMK 
GOTO  112 

ELSEIF  (LINE.EQ.7)  THEN 

PRINT  x.'EtITER  THE  EMISSIVITY  OF  THE  LID.’ 

READ  X,  EPSL 
GOTO  112 

ELSEIF  (LINE.EQ.A)  THEN 

PRINT  X, 'ENTER  THE  EMISSIVITY  OF  THE  DIE.’ 

READ  X,  EPSD 
GOTO  112 

ELSEIF  (LINE.EQ.9)  THEN 

PRINT  x.'DO  YOU  WISH  TO  CHANGE  YOUR  CHOICE  OF  HEAT  INPUT?’ 

READ  (X.A)  ANS 

FORMAT  (All 

IF  (ANS. EQ. YES)  THEN 

GOTO  110 

ELSE 

PRINT  X, 'ENTER  NODE  AND  HEAT  INPUT  TO  BE  CHANGED!* 

HEAD  X.  CMNODE,  CHIIATT 
IKGMMtlDE,  //())=99')1 
A(CMMODE,  7  70  )=(:MWATT 
PRINT  x.'ALL  CHANGES  COMPLETED?’ 

READ  (X,6)  CHANS 
FORMAT  (Al) 

IE  (CHANS . EQ .YES)  THEN 
CHOICE =3.0 
GOTO  112 
ELSE 

GOTO  114 
ENDIF 
ENDIF 

ELSEIF  (LINE.NE.  1  . AND . L I  ME . NE . 2 . AND . LI  ME . NE , 3 . AND , L INE . ME . 4 . AND . 
•FLINE.NE.5. AND. LINE.NE.6 .AND.LINE.NE.7 .AND.LINE.NE.S .AND. 
+LINE.NE.9)  THEN 

PRINT  X, ' try  AGAIN’ 

GOTO  112 
ELSE 

CONTINUE 

ENDIF 


Kl  =  1101  . 
K2=1102. 
K3=1103. 
K4=110A. 
K5=ll()5. 
K6  =  lin6  . 
K7  =  1107  . 
K8=1108. 


K1  -  K8  ASSIGN  THE  CURVES  llOl  -  1108  TO  SPECIFIC  NODES  FOR  USE 
IN  THE  THERMAL  ANALYZER. 


THE  FOLLOWING  SETS  OF  ALGORITHMS  ESTABLISH  THE  NODE  AND  THE 
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MODES  HHICM  ARE  COdtlECTED  TO  IT.  THEY  ALSO  ESTABLISH  THE  CURVE  TO  BE 
USED  IN  HIE  ANALYSIS. 


113  DO  31  Y  =  l,YROll 

DO  32  X=l,XRnU 
C 

I=COUIMER 

C 

IF  (X-1 .EQ.O.AMD.Y-1 .EQ.O)  THEM 
A(I,I*n=  K1 
A(I.I+HR)=  K2 
B(I,I+1)=  10»(I+1)+1 
D(I,I+MR)=  10m(I+HR)+1 
C 

ELSEIF  (Xtl .GT.XROW. AMD, Y-I. EQ.O)  THEM 
A(  I,  I-l  )=  K1 
A(I,I+MR)=  K2 
B(  I .  I-l)=  10K(  I-n  +  1 
D(1.I+MR)=  10k(I+MR)+1 
C 

ELSEIF  (Y-l.EQ.O)  THEN 
A(I,.I-1)=  K1 
A{  I ,  I  +  l  )=  KI 
A(I,I+MR)=  K2 
D(I,I-l)=  10*(l-l)+l 
B(I,I+1)=  10X(I+l)+l 
B(I,I  +  tlR)=  10x(I+MR)  +  l 
C 

ELSEIF  (X-1. EQ.O. AMD. Y+l.GT.YROlO  THEM 
A(I,I+1)=  K1 
A(l,I-MR)=  K2 
B(I,IU)=  10x(I  +  n  +  l 
n(I,I-MR)=  lOXd-MIO  +  l 
C 

ELSEIF  (X+1 .GT .XROH.AMD.Y+l .GT.YROU)  THEM 
A  (  I ,  I  - 1  )  =  K 1 
A(I,I-MR)=  K2 
B(  I ,  i-l)=  10x(  I-n  +  1 
D(I,I-MR)=  10x(I-MR)+l 
C 

ELSEIF  (X-1. EQ.O)  THEM 
A(I,I+1)=  Kl 
A(I,I-MR)=  K2 
A(I,I+MR)=  K2 
D(  I,  I  +  l  )=  10x(  I  +  n  +  I 
0(1,1-110)=  10X(I-HR)  +  1 
B(I,I+MR)=  10X(I+MR)+1 
C 

ELSEIF  (X+1 .GT.XROH)  THEM 
A(  I, I-l  )=  Kl 
A(I,I-MR)=  K2 
A(I,I+MR)=  K2 
B(I,I-1)=  10X(I-1)+1 
B(I,I-MR)=  10X(I-HR)+1 
B(I,I+MR)=  10X(I+MR)+1 
C 

ELSEIF  (Y+1 .GT.YROU)  THEM 
A(I,I-1)=  Kl 
A(I,I+1)=  Kl 
A(I,I-MR)=  K2 
B( I ,  I-l  )=  10x( I-l )  +  l 
B(I,I+1)=  10x(I+l)+l 
B(I,I-MR)=  10X(I-MR)+1 
C 

ELSE 

AC  I ,  I-l  )=  Kl 
A( I , I  +  l  )=  Kl 
A(I,I-MR)=  K2 
A(I,I+MR)=  K2 
B(I,I-1)=  10x(I-l)+l 
D(I,I+1)=  10x(I+l)+l 
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c 

52 

51 


V(  1  .  l-IIKI* 
S(  I.  l*IIR)t 


£IIDIF 


KS 

Ds  I .  I  Mil  )  >  1  Ui(  1  Mil  )«  I 
COUtirER'COUMTERf  I 


coiiTitmc 

COtITlilUE 


KHTRMllM  '  ' 

DO  55  Ysl.YROH 
DO  54  X*l.XROH 

1=KHTR 

IF  (X-l.EQ.O.AHD.Y-l.EQ.O)  THEM 
A<  1 .  1*1  )s  X4 
A(  I.  IMIR):  K5 
B(I.1«1)=  10i(t4n«l 
B(l.l4ltR)=  10)l(l4t|R)4l 

ELSEIF  (X>1.EQ.O.AMD.Y41.6T.YROII>  THEM 
A(  1.14  1):  K4 
A(  I .  l-IIR):  KS 
8(1.141):  104(141)41 
B(l.l-MR):  10H(1-I|R)41 

ELSEIF  (X41 .cr.XROH.AMD.T^l.EQ.O)  THEM 
Ad. 1-1):  K4 
Ad.lMIR):  KS 
R(I.I-1>>  104(1-1)41 
Bd.lMIR):  10l(l4MR)4l 

ELSEIF  (X4l.G1.XR0M.AIID.Y4l.CT.YR0H)  THEM 
.  A( 1,1-1):  K4 

A((.I-MR):  KS  ..S' 

Bll,I-l)«  10«(I-l)4l  .;V, 

B(I,1-MR):  10K(l-MR)4l 

ELSEIF  (X-l.EO.O)  THEM 
A( 1 . 14l  )>  k4 
A(I.I-IIR)s  KS 
A(1.I4||R):  KS 
Bd,  141)3  101(141)41 
Bd.l-MR):  10i((-IIR)4l 
Bd.dllK):  10l(l4MR)4l 

ELSEIF  (X41.GT.XR0M)  THEM 
Ad, 1-1):  K4 
Ad.I-MR):  KS 
A(  I.  I4||R)3  KS 
B(I,I-1):  104(1-1)41 
Bd.I-MR):  10K(I-MR)4l 
B(I,I4||R)X  10l(l4|IR)4l 

ELSEIF  (Y-l.EO.O)  THEM 
Ad  ,  1-1  )•  K4 
A( 1 . I4l  )3  K4 
A(1,I4||R)*  KS 
Bd.I-l):  104(1-1)41 
Bd. 141)1  101(141)41 
Bd.lMiR):  104(IMIR)4l 

ELSEIF  (Y4l.CT.YR0M)  THEM 
Ad.I-l):  K4 
Ad. 141)3  K4 
Ad.l-IIR):  KS 
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D(i,i-i)=  ioK(r-i)+i 

n  ( I ,  I  f  1 )  =  1 0  N  ( I  + 1 ) » 1 

R(I,I-IIR)=  10k(I-IIR)  +  1 

c 

ELSE 

A(  I ,  I-l  )= 

A(  I ,  I  n  )=  K‘'( 

A(1.I-HR)=  K5 
A(I,I  +  »IR)=  K5 
D(I,I-1)=  10K(I-1)+1 
n(  I ,  I  +  l  )=  10X(  I  +  l  )♦  I 
l)(l,I-tlR)=  10K<I-||R)  +  1 
D(I.I+HR)=  10X(I+UR)+1 
ENDIF 
C 

A(1,I-?JL)=  K3 
B(I,I-ML)=10*(I-NL)+1 
C 

KNTR  =  KflTR+l 
C 

34  COMTiriUE 

33  COUTIIIUE 

C 
C 

CflTR  =  2»<ML  +  l 
HP  =  IIL  +  1 

DO  35  Y  =  1,IIY 
DO  36  X=1,MX 
C 

I=CMTR 

C 

IF  (X-1  .EO.O.AIlD.Y-1  .EQ.O)  THEH 
A(  I,  I  +  n=  K4 
A(  I  ,  I  i(IX)=  K5 
11  (  1  ,  1  t  1  )  =  1  0  >  ( I  ♦  1)  t  I 
U(I,HNX)=  10x(I  +  l)X)  +  l 

c 

ELSEIF  (X-1 .EQ.O. AND, Y+l.GT.HY)  THEM 

A(I.l+n=  KA 

A(  I  ,  I-|IX)=  K5 

11  ( I .  I  n )  =  1 0  K  <  1 1  n  n 

D(1,1-MX)=  10X(I-MX)+1 
C 

ELSEIF  (X+1 .GT.MX.AND.Y-1 .EQ.O)  THEN 
A(I,I-1)=  K4 
A(I,I+NX)=  K5 
0(  I  ,  I-l  )=  10K( I-l )  +  l 
D(I,I  +  NX)=  lOxd  +  NXJ  +  l 
C 

ELSEIF  (X+1 .GT.NX.AND.Y+1 .GT.NY)  THEN 
A(I,I-1)=  K4 
A(  I  ,  I-NX)=  K5 
B(  I,  I-l  )=  10X(I-1  )  +  l 
D(I,I-NX)=  10*(I-NX)+1 
C 

ELSEIF  (X-l.EQ.O)  THEN 
A( I , I  +  l  )=  K4 
A(I.I-NX)=  K5 
A(I,I+NX)=  K5 
B(I,I+1)=  10X(I+1)+1 
B(I,I-NX)=  10«(I-MX)+1 
D(I,I»NX)=  10»c(I  +  NX)  +  l 
C 

ELSEIF  (X+l.GT.NX)  THEN 
A(  I, I-l)=  K4 
A(I,I-NX)=  K5 
A(  I  .  I  •NX)=  K5 
U{  I  ,  I-l  )=  10»«(  I-l  )  +  l 
D(I,I-NX)=  10X(I-NX)+l 
B(I,I+NX)=  10X(I+NX)+1 
C 

ELSEIF  (Y-1 .EQ.O)  THEN 
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A(I,I-n=  K4 
AC  I.  HI  )=  K^t 
A(l,l+nx)=  K5 

D(i,i  +  n=  io*(i  +  i)  +  i 

B(I,1  +  UX)=  10*(1  +  I1X)  +  1 

ELSEIF  (Ytl.GT.HY)  THEH 
A(I,I-l)=  KA 
A(I,I  +  1)=  K-i 
A(I,I-riX)=  K5 
B(I,I-1)=  10K(I-1)+1 
BC1,I  +  1)=  loxci  +  n  +  i 
B(i,i-nx)=  io*(i-nx)+i 

ELSE 

A(I,I-1  )=  K4 
A(  I .  I+l  )=  K'i 
A(I.I-IIX)=  K5 
A(I.I  +  IIX)=  K5 
B(  I.  I-l  )=  )  +  l 

B(l,l  +  1)=  10X(I  +  1  )♦! 
B(I,I-tlX)=  10>c(I-riX)  +  l 
B(I,I  +  IIX)=  10K(I  +  IJX)-H 
ENDIF 

ALPHA  =  -IIP+X+(Y-n*(tlR  +  riX) 

BETA  =  -(IL  +  X  +  (Y-1  )»(I(R  +  HX) 

GAHHA  =  -tlPUKY-l  )«( IIR*IIX)4||R 
DELTA  =  -IIL  +  X+(Y-1  )x(IIR  +  NX)  +  HR 

A(I,I+ALPHA)=  K6 

AC  I,  I+BETA)=  K6 

AC  I  ,  HCiAM(1A)s  K6 

AC  I , I ' DEL  I  A )=  K6 

AC  1  ,  1  Mil  /A  )=  K7 

AC  I  tALPHA,  n=  K6 

ACHI)ETA,I)=  K6 

AC  I  tCAIlIlA,  n=  K6 

AC I  +  DELTA, I  )=  K6 

DC  I , I  +  ALPHA)=  10»( I  +  ALPHA)  +  1 

DCI,nDETA)=  10K<  I  +  BETA)  +  1 

DC  I ,  HGAMI1A)=  10K(  l+GAtniA)  +  l 

DC  I. I  +  DCLTA)=  10»( I  +  DELTA)  +  1 

DC  I,  I  +  (IL/A)=  10«(  I+tlL/A)  +  l 

DC  14  ALPHA, I  )=  10»I  +  1 

DCI+DETA,I)=  lOXl+1 

DC  HGAHHA,  I  )=  10*1  +  1 

B(1+DELTA,I)=  10*1+1 

CHTR  =  CIITR+1 

CONTimiE 

COMTIHUE 


CTR=2*HL+ML/A+1 
DO  37  Y=1,MY 
DO  38  X=1,NX 

I=CTR 

IF  (X-l.EQ.O.AHD.Y-l.EQ.O)  THEM 
ACI,I+1)=  KA 
ACI,I+NX)=  K5 
DCI.I+1)=  10*CI+1)+1 
UCI,l4tlX)=  10*<I  +  HX)  +  1 

ELSEIF  (X-1 .EQ.O.AMD.Y+l.GT.MY)  THEN 
ACI,I+1)=  KA 
ACI,I-NX)=  K5 
DCI,I+1)=  10x(I+l)+l 


45 


B(I,I-UX)=  10K(I-HX)+I 

ELSEIF  (X+1 .GT .nx.AND.Y-1 .EQ.O)  THEN 
A(I,I-1)=  K4 
A(I,1+NX)=  K5 
D(I,1-1)=  +  1 

n(i,iuix)=  )o<<M  +  nx)*i 

ELSEIF  (X+1  .GT.IJX.AND.Y+1  .GT.liY)  THEM 
A(r.I-l)=  K<« 

A(I,I-IIX)=  K5 
D(I,I-1)=  10X(I-1)+1 
B(I,I-tlX)=  10K(I-UX)  +  1 

ELSEIF  (X-l.EQ.O)  THEN 
A(  I,  I  +  l  )=  K<i 
A(I,I-NX)=  K5 
A(I,l  +  tlX)=  K5 
B(I,1+1)=  10X(I+1)+1 
D(I,I-NX)=  10K(1-NX)+1 
D(1,1+NX)=  10*(I+NX)+1 

ELSEIF  (X+l.GT.MX)  THEN 
A(I,I-1)=  KA 
A(I,I-NX)=  K5 
A(I,I  +  tlX)=  K5 
B(I,I-1)=  10x(l-n  +  l 
B(I,I-NX)=  10X(I-NX)+1 
B(I,I+NX)=  10X(I+NX)+1 

ELSEIF  (Y-l.EQ.O)  THEN 
A(  I,  I-l  )=  K<* 

A(  I  ,  I*  n=  Ky. 

A(  1  ,  I  MIX)=  K5 

U<t.I-n=  1«>.<I-1)  +  1 

1)(I,I  +  1)=  1UX(I  +  1)  +  1 
B(I,I+NX)=  10x(I+NX)+l 

ELSEIF  (YU.GT.NY)  THEN 
A( I , I -1  )=  KA 
A(  I , I  +  l  )=  KA 
A( I , I-NX)=  K5 
D(I,I-1)=  10X(I-1)+1 
B(I,I+1)=  10X(I+1)+1 
D(I,I-NX)=  10X(I-NX)+1 

ELSE 

A(I,I-n=  KA 
A( I , I  +  l  )=  KA 
A(I,I-NX)=  K5 
A(I,I+MX)=  K5 
B(I,I-1)=  10X(I-1)+1 
D(I,I+1)=  10X(l+l)+l 
B(I,I-NX)=  10X(I-NX)+1 
•  B( I , I+NX)=  10X( I+NX)+1 

ENDIF 

A(I,I-NL/A)=  K7 
B(I,I-NL/A)=  10x(I-NL/A)+1 

CTR  =  CTRH 

CONTINUE 

CONTINUE 


TH2=(  TIIICK-THl  )/(LVL-l  ) 
DELX  =  l  EUGril/XROU 
DELYMIinTH/YROll 
DELZ1  =  TII1 
DELZ2=TH2 
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DELZ3=(THl+TH2)/2 
DEL2<i  =  Ttl2/2 
KAIR=1109. 

C 

C 

C  THE  FOlintlTrir,  ALOORITHM  hill  calculate  a  CORRECTIOri  FACTOR  FOR  THE 
C  THERMAL  COIIDUCT  1 V  i  1  f  01  AIR  TO  INCLUDE  RADIAIlOtl  CIFECIS.  THIS  IS 
C  DOME  TO  SAVE  OH  COMPUTER  RUM  TIME  IM  THE  THERMAL  ANALYZER.  SEE 
C  THESIS  FOR  GREATER  DETAILS. 

C 

C 

IF  (IlNIT.EQ.l)  THEN 
LIUK=LIU+273.15 
ELSE 

LIDK=(5/9)x(LID-32)+273.15 

ENDIF 

C 

EMS=(l/EPSL)+(l/EPSD)-l 

FE=1/EH5 

TX=250.0 

DO  200  N=l,9 

TAV=0.5X(LIDK+TX) 

THAIR{N)=7 .76E-5XTAV+2.87E-3 
HGAP=THAIR(N)/GAP 

HRAI)=5.669E-8kFEx(LIDK»»2+TXkk2)«<LIDK+TX) 

RAT=HRAD/HGAP 

FACT0R=1+RAT 

KAR(II)=THAIR(N)XDELXKDELYKFACT0R*0. 0001/GAP 
TX=IX+50 
200  CONTINUE 


THE  roiinillNfi  GROUP  OF  ALGORITHMS  CONNECT  THE  TOP  LAYER  NODES  TO 
THE  LID  MODES  (NODES  751-765). 

THE  THERMAL  ANALYZER  IS  LIMITED  TO  A  MAXIMUM  OF  27  CONNECTIONS. 
THEREFORE  DY  CONNECTING  THE  TOP  LAYER  NODES  IN  GROUPS  OF  TWENTY  THIS 
LIMITATION  CAN  BE  BYPASSED.  THE  SAME  ANALOGY  IS  APPLIED  TO 
CONNECTING  THE  BOTTOM  LAYER  NODES  TO  THE  SINK  NODES  (NODES  766-769). 


62 

C 


63 

C 


64 

C 


65 

C 


66 

C 


67 

C 


DO  62  1=1,20 

A( I,75l)=KAIR 
B(I, 7511=7511 
CONTINUE 

DO  63  1=21,40 
A(I,752)=KAIR 
D(I, 7521=7521 
CONTINUE 

DO  64  1=41,60 
A(  I,753)  =  KAIR 
B(I,  7531=7531 
CONTINUE 

DO  65  1=61,80 

A(  I,754)=KAIR 
B(I,  7541=7541 
CONTINUE 

DO  66  1=81,100 
A(  I  ,  7  55) =KAIR 
B(  1 ,7551=7551 
CONTINUE 

DO  67  1=101,120 
A(  I  ,  756  )=KAIR 
D(  1 ,  756  1=7561 
CONTINUE 

DO  68  1=121,140 
A(I,757)=KA1R 
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D(I,757)=7571 

68  COMTltlUE 
C 

DO  69  I  =  160 

A( I ,758 )=KA1R 
D(I.753)=7581 

69  CONIIlilir 
C 

DO  70  1=161,180 
A(1,759)=KA1R 
B(l,759)=7591 

70  CONTINUE 
C 

DO  71  1=181,200 
A( I,760)=KAIR 
B(  I,  760)  =  7601 

71  CONTINUE 
C 

DMri=2xNL  +  NL/A+l 
C 

DO  85  1=DMH,DMM*20 
A(  1 ,766  )  =  K7 
B(I,766)=7661 
83  CONTINUE 

C 

DO  89  I-DMM+21, DMM+90 
A(  1 , 767  )  =  K7 
B(I,  767  1=7671 
89  CONTINUE 

C 

IF  (AN51I.EQ.1.AND.ANSII.EQ.9)  THEN 
LI1)N=13 
SNKN=9 

DO  72  1=201,220 
A(  I  ,761  )  =  KAIR 
IK  1  .  761  )  =  761  1 

72  CONTINUE 
C 

DO  73  1=221,290 
A(  I  ,762  )=KMR 
U(  1  ,  762  1  =  7621 

73  CONTINUE 
C 

DO  79  1=291,256 
A(  I  ,763)=KAIR 
IK  1,763  1  =  7651 
79  CONTINUE 

C 

DO  8  5  I  =  DI1M+91 ,  DMM+6  0 
AC  1 ,7681=K7 
6(1,7681=7681 

85  CONTINUE 
C 

DO  86  I=DHM+61,DH1U69 
A(I,7691=K7 
B(  1,7691=7691 

86  CONTINUE 
C 

ELSEIF  (ANSII.EQ.2)  THEN 
L I DN=15 
5NKN=9 

DO  75  1=201,220 
AC  1 , 76 1 1=KAIR 
DC  1,761  1  =  761  1 

75  CONTINUE 
C 

DO  76  1=221,290 
AC  1 ,7621=KAIR 
0(1,7621=7621 

76  CONTINUE 
C 

DO  92  1=291,260 
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92 

C 


95 

C 


94 

C 


87 

C 


88 

C 


77 

C 


78 

C 


79 

C 


80 

C 


81 

C 


89 

C 


90 

C 


82 

C 


A( I , 765)=KAIR 
IK  ] ,  765)  =  7631 

coin  imiE 

DO  93  1=261,280 
A( I , 764)=KAIR 
IK  I  ,  764)  =  76'*1 
com  lllUL 

DO  94  1=281,288 
A(  I  ,  765)=KAIR 
D(  1  ,  765)  =  7651 
CONTI  HUE 

DO  87  I=DMH+41 , DMMf60 
A(I,768)=K7 
B(  1  ,  768  )=7681 
CONTINUE 

DO  88  I  =  DIUn61 ,  DMH+72 
A(  1 , 769  )  =K7 
B(I,769)=7691 
CONTINUE 

ELSEIF  (AN5II.Eq.3)  THEN 
LIDN=15 
SNKN=4 

DO  77  1=201,220 
Ad, 761  )=KAIR 
IK  I , 76  1  )  =  76  1  1 
CONTINUE 

DO  78  1=221,240 
At  I,  762  )=KAIR 
IK  I  ,  762  )  =  7621 
CONI  IIIUE 

DO  79  1=241,260 
A( I,763)=KAIR 
B< I ,7631=7631 
com  INUE 

DO  80  1=261,280 
A(I,764)=KAIR 
D( 1,7641=7641 
CONTINUE 

DO  81  1=281 , 300 
Ad  ,  765)=KAIR 
Bd,  7651=7651 
CONTINUE 

DO  8  9  I  =  DMrn-41,DMM+6  0 
■Ad,7681  =  K7 
BCI, 7681=7681 
CONTINUE 

DO  90  I=DMM+61 , DMM+75 
Ad,769  1=K7 
Bd  .  769  1=7691 
CONTINUE 

ELSEIF  (AN51T.CQ.51  THEN 
LIDN=ll 
SNKN=3 

DO  82  1=200,216 
Ad  ,  76  I  )=KAIR 
IK  I  ,761  1  =76  1  1 
CONTINUE 

DO  91  I  =  DHHf41,DMri+54 
A( I ,768  1=K7 
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ooooooooo  oo 


D( I ,768)=7681 
91  COHHNIIE 
C 

EtlDIF 

C 

C 

DO  16  I  - 1 ,  L  I  0(1 
L  I  III  (  I)=L  ID 
16  CONTIIIUE 

C 

DO  19  I  =  l.SriK(l 
SIIIKT(I)=5iriK 
19  COIITIMUE 

c 

c 

C  INITIAL  SETTINGS  FOR  THE  THERMAL  ANALYZER 

C 

C 

HTR  =  0 
EXP  =  0 
SECIir  =  0 
C0EFCV=9 
FASTIIT  =  0 
HTCV  =  0 
UNK  =  0 

CONTROL1=0 

C0NTR0L2=0 

CONTROL3=0 

NnDC5=300 

KMAX=50 

HMAX=6 

XP  =  1 

SllT-2 

El  CV  - 1  5 

sun  ^:> 

TCV  =  8 

ACC= . 05000 
DAMP= .6667 
MAXITS=12 
C()llvrAC=  .8000 


C(1)=(DELYXDELZ1 )/DELX 
C(2)  =  (DELXXDELZ1  )/DELY 
C( 3)=( DELXKDELY)/DELZ3 
C(A)=(DELYKnELZ2)/DELX 
C(  5)  =  (  I)ELXXDELZ2)/DELY 
C(6  )  =  ( DELXXDELYl/DELZZ 
C(7  )  =  ( AXDELX*DELY)/DELZ2 
C(8)=( DELXXDELYI/DELZA 


THE  TEMPERATURES  SELECTED  FROM  THE  SILICON  THERMAL  CONDUCTIVITY 
CURVE.  THE  FIRST  GROUP  IS  IN  DEG  CENTIGRADE,  THE  SECOND  IN 
FAHRENHEIT.  EACH  TEMPERATURE  GROUP  IS  FOLLOWED  BY  THE  THERMAL 
CONDUCTIVITY  OF  SILICON  AT  THAT  TEMPERATURE  IN  EITHER  WATTS/CM-C 
OR  DTU/HR-IN-F. 


IF  (UNIT.EQ.l)  THEN 
TEMPI =-271 .2 
TEIU’2  =  -267 . 2 
TtMI'3  =  -263.2 
TEMP6=-255.2 
TEMP5=-2A3.2 
TEfir6=-228 ,2 
TCMI’/  =-208 . 2 
T  Kill’s  =  -17  3. 2 
TEMP9=26 .85 
THCOMIK  1 )=0 .25 
THCOIII)(2)=6.5 
THCUIII)(3)  =  10.0 
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THC0tlD(<i)  =  17 .0 
7HCOtin(5)  =  19.0 
THCOIII)(6)  =  17.0 
THCOIin(7)  =  10.0 
TMC0t)D(8)  =  5.5 
THC0IIU(9)  =  1 .5486 
Eir.l 

1I;MC  I  -4  56.(1/ 

Tl:MI*2  =  -448.87 
TENP3=-44l .67 
TEMP4=-427 .27 
TEMP5=-405.67 
TEHP6=-378 .67 
TEMP7=-342.67 
TEMP8=-279 .67 
TEHP9=80.33 
THCOtlD(l)  =  1.20 
THC0I1D(2)  =  21 .67 
THC0UD(3)=48.15 
THC0MD(4)=81 .85 
THC0HD(5)=91 .48 
THC0IID(6)=81 .85 
THCOIID(7)=48.15 
THC0lin(8)  =  26 .48 
THC011D(9)=7 .46 
EIIDIK 
C 
C 

DO  23  1=1,8 
nn  24  J=l,9 

(:URVE(  I,  J)=C(1  )KTHCOHD(  J) 

24  COfITKIUE 

23  COtITIIIUE 

C 

DO  39  1=1,9 
5ll()l(I)=0 
SUDC I  )=0 
39  COIITIIIUE 

C 

c 

C  OPEfI  A  FILE  FOR  THE  OUTPUT.  THIS  FILE  HILL  BE  111  THE  CORRECT 

C  FORMAT  TO  BE  READ  BY  THE  THERMAL  AHALYZER. 

C 

C 

OPEN  ( 3, FI LE='NICCKT' ,FORM=* FORMATTED* ,ACCESS=* DIRECT ' , 

+  RECL  =  108,5]ATU5='NEII' ) 

C 

HRITE  (3,701)  TITLE 

701  FORMAT  <1X,A79) 

WRITE  (3,702)  R,  CONST,  HTR,  EXP,  SECHT,  COEFCV,  FASTHT,  HTCV, 
+UIIK,  UNIT 

702  FORMAT  (1316) 

WRITE  (3,703)  CONTROLl,  C0IITR0L2,  CONTROLS 

703  FORMAT  (316) 

WRITE  (3,704)  NODES,  KMAX,  HHAX,  XP,  SHT,  EFCV,  STHT,  TCV,  UNK, 
+  UNIT 

704  FORMAT  (1316) 

WRITE  (3,705)  ACC,  DAMP,  MAXITS,  CONVFAC,  INTEMP 

705  FORMAT  (2GI2. 4, 112,2012.4) 

C 

DO  25  1=1,8 

Winil-  (3,7(16)  TEMPI,  CURVEd.n,  TEMP2,  CURVE(I,2), 

+  1  I  MP .T, CURVE!  I,  3),  TEMP4,  CURVE(I,4),  TEMPS,  CURVE!  1,5),  TEMP6, 
+CURVE(I,6),  TEMP7,  CURVE(I,7),  TEMP8,  CURVE(I,8),  TEMP9, 

+CURVE( 1,9) 

706  FORMAT  (9012. 4) 

25  COliriNOE 

C 

DO  201  N=l,9 

WRITE  (3,710)  THAIR(N),  KAR(H) 

710  FORMAT  (9012.4) 

201  CONTINUE 
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HRITE  (3,707  )  (LIDT(H),H=1,12),  ( SIHKT(H) , N=1 , 2 ) 

707  FORMAT  (9012.4) 

C 

DO  50  1  =  1, R 

COUtlT  =  0 
COUIirER  =  0 
DO  fii  j  =  i.aon 

If  (A(I.J).EQ.O)  THEN 
GOTO  51 
ELSE 

C0UHT=C0UNT+1 
TEHP(COUriT)=A(I,  J) 

END  IF 

51  COflTinOC 

DO  52  J=l,800 

IF  (B(I,J).EQ.O)  THEM 
GOTO  52 
ELSE 

COUMTER=COUHTER+l 
THP(COUIITER)  =  B(l,  J) 

EHDIF 

52  CONTINUE 
C 

WRITE  (3,708)  COUNT,  ( THP( N) , H=1 , COUNTER) ,  ( SUBl ( N ) , N=COUNT+l , 9 ) 

708  FORMAT  (15,17,8112) 

WRITE  (3,709)  ( TEMP( N ), N=1 , COUNT) ,  (SUB( N) , N=COUNT4l , 9 ) 

709  FORMAT  (9G12.4) 

50  CONTINUE 

C 

C 

ENDFILE  (3) 

CLOSE  (3) 

C 

C 

STOP 

END 
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APPENDIX  B 


The  following  will  demonstrate  a  typical  run  using  MICCKT.  On  screen  prompts  will 
be  indicated  within  brackets  (  <  >  ),  responses  (Resp)  will  be  directly  afterwards. 
Comments  will  be  inserted  for  clarity. 

<  Enter  the  name  of  your  file  > 

Resp:  This  is  an  alphanumeric  string  with  a  maximum  of  76  characters. 

<  You  may  select  one  of  the  following  options; 

1. )  16  by  16,  yields  640  nodes 

2. )  12  by  24,  yields  720  nodes 

3. )  10  by  30,  yieds  750  nodes 

4. )  8  by  32,  yields  640  nodes 

5. )  12  by  18,  yields  540  nodes  > 

Resp:  3 

<  Enter  1.)  for  SI  or  2.)  for  English. 

Resp:  1 

<  Input  chip  dimensions  in  centimeters. 

Enter  length:  > 

Resp;  0.75 

If  the  English  system  of  measurement  was  chosen  the  dimension  prompt  would  be  in 
inches. 
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<  Enter  width:  > 

Resp:  0.5 

<  Enter  thickness:  > 

Resp:  0.05 

This  will  be  the  overall  thickness  of  the  chip.  The  top  layer  is  automatically  set  at 
0.002  inches  or  0.00508  centimeters.  Subsequent  layer  thicknesses  are  derived  as  follows: 
(overall  thickness  -  first  layer  thickness)/3. 

<  Enter  the  thickness  of  the  air  layer  between  the  lid  and  the  top  of  the  die.> 

Resp:  0.2 

The  maximum  allowable  thickness  is  0.254  cm  (0.10  in). 

<  Enter  the  temperature  of  the  lid:  > 

Resp:  40 

Temperatures  will  be  in  Centigrade  or  Fahrenheit  depending  on  the  users  choice  of 
dimensions. 

<  Enter  the  temperature  of  the  sink:  > 

Resp:  15 

<  Enter  the  emissivity  of  the  lid:> 

Resp:  0.7 

<Enter  the  emissivity  of  the  chip:> 

Resp:  0.8 
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<  Heat  inputs  will  be  in  watts.> 

This  prompt  is  a  reminder  of  the  units  the  user  is  working  in.  Heat  inputs  will  be  in 
Btu/hr  if  in  the  English  system. 

<  Heat  input  to  the  microcircuit  can  be  entered  in  the 
following  manner.  Choose  1.)  to  enter  as  total  heat 
applied  to  the  chip,  2.)  to  enter  as  the  average  heat 
per  area,  3.)  to  enter  node  by  node,  or  4.)  for  no  heat 
input. 

Resp:  2 

<  Enter  the  average  heat  per  unit  area:> 

This  will  be  in  either  Watts/cm^  or  Btu/hr-in^. 

Resp:  5 

<  Summary  of  input  data. 


1.)  Length  of  the  chip  is  (cm): 

0.7500 

2.)  Width  of  the  chip  is  (cm): 

0.5000 

3.)  Thickness  of  the  chip  is  (cm): 

0.5000E-01 

4.)  Thickness  of  the  air  layer  (cm): 

0.2000 

5.)  Temperature  of  the  lid  is  (deg  C): 

40.00 

6.)  Temperature  of  the  sink  is  (deg  C): 

15.00 

7.)  The  emissivity  of  the  lid  is: 

0.7000 

8.)  The  emissivity  of  the  chip  is: 

0.8000 
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0.6250E-04  Watts  > 


9.)  The  heat  inputs  to  the  node  are: 

Heat  applied  to  each  node  is: 

Heat  applied  was  5  Watts/cm^,  this  was  multiplied  by  the  area  of  the  chip  and  then 
divided  by  the  number  of  nodes  in  the  top  layer  to  produce  an  average  heat  input  to  each 
node.  If  option  1  had  been  chosen  the  total  heat  applied  to  the  chip  would  be  divided  by  the 
number  of  nodes  in  the  top  layer  to  get  an  average  input  to  each  node. 

<  Do  you  wish  to  make  any  changes,  y  or  n?> 

Resp:  y 

<  Which  line  item  do  you  wish  to  change?> 

Resp:  4 

<  Enter  the  thickness  of  the  air  layer  between  the  lid  and  the  top  of  the  die.> 

Resp:  0..25 

<  Summary  of  input  data. 

1. )  Length  of  the  chip  is  (cm): 

2. )  Width  of  the  chip  is  (cm): 

3. )  Thickness  of  the  chip  is  (cm): 

4. )  Thickness  of  the  air  layer  (cm): 

5. )  Temperature  of  the  lid  is  (deg  C): 

6. )  Temperature  of  the  sink  is  (deg  C): 

7. )  The  emissivity  of  the  lid  is: 

8. )  The  emissivity  of  the  chip  is: 


0.7500 

0.5000 

0.5000E-01 

0.2500 

40.00 

15.00 

0.7000 

0.8000 
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9.)  The  heat  inputs  to  the  node  are: 

Heat  applied  to  each  node  is:  0.6250E-04  Watts  > 

Do  you  wish  to  make  any  changes,  y  or  n?> 

Resp:  n 

With  a  no  response  the  program  will  generate  an  output  file  for  use  with  the  thermal 
analyzer. 
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